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Optical signal processing can be performed in either digital 
or analog systems. This dissertation is concerned with the use of 
digital computers and coherent optical systems in optical signal proc- 
essing, such as image restoration, phase-object visualization, image 
contrast reversal, optical computation, image multiplexing, and 
fabrication of spatial filters. 

Part I of the dissertation is concerned with digital optical data 
processing. It deals with restoration of images degraded by 
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signal-dependent noise. When the input data of an image restoration 
system are the numbers of photoelectrons received from various 
areas of a photosensitive surface, the data are Poisson distributed with 
mean values proportional to the illuminance of the incoherently radiat- 
ing object and background light. In this model, the shot noise is predomi 
nant in causing random variation of the observed data. The noise is 
signal-dependent. A nonlinear restoration scheme based on the prin- 
ciple of maximum likelihood is derived, which is optimum for a broad 
class of objects with preassumed average contrast and size of typical 
details. The nonlinear estimation scheme can be solved by linear 
iteration, or methods of nonlinear programming. 

One -dimensional objects distorted by diffraction or motion blur 
are estimated in computer simulations. The performance of the esti- 
mator is measured by percentage squared error. When the background 
illuminance is low, the maximum- likelihood estimate is much 
superior to the linear least-squares estimate. 

Part II of the dissertation is concerned with optical signal 
processing using coherent optical systems. Following a brief review 
of the pertinent details of Ronchi's diffraction grating interferometer, 
moire effect, carrier -frequency photography, and achromatic holog- 
raphy, two new shearing interferometers based on them are presented. 
Both interferometers can produce variable shear. 
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New methods for optical testing and phase-object visualization, 
spatial differentiation, complex amplitude addition and subtraction, 
image contrast reversal, two-dimensional photometry, image multi- 
plexing, and fabrication of spatial filters can be achieved easily by 
using these two new interferometers. All the new methods are com- 
pared with the existing methods. They all have some definite advan- 
tages in either flexibility, simplicity, light efficiency, or real-time 
operation. Finally, a new method for focusing and alignment is also 
developed. 
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Introduction 


Because of the two-dimensional nature of pictures, an optical 
communication system can transmit much more information than 
electronic systems per unit time interval. An optical system can be 
thought of as a bank of parallel electronic communication systems. Due 
to the problem of the already over-crowded radio communication chan- 
nels, much effort has been put into the study of optical communication 
systems, which afford wide temporal bandwidth and multiple spatial 
channel capacity. Optical signal processing is as important in optical 
communications as its counterpart in electronic communications. In 
this thesis we will discuss methods of optical signal processing using 
both digital computer and coherent optical systems. 

In part I of this thesis we will present methods for the restora- 
tion of images degraded by signal- dependent noise. Processing of radio 
signals (deterministic or stochastic) degraded by a radio communica- 
tion system is based on a linear model of signal formation, in which 
the noise is independent of the desired signal and combined with it 
additively; the sum of these constitutes the observed data. In signal 
estimation, the data are passed through a linear spatial filter whose 
transfer function is designed to reproduce the original signal within the 

smallest possible mean-square error. This kind of filter is known as 

( 1 - 3 ) 

the Wiener filter or matched filter. This additive -noise model 

has been used extensively with success. 
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The signal processing techniques used at radio frequencies have 
been extended to the optical frequency region. Methods for restoring 
images degraded by an optical system are often based on a linear model 
of the optical process: a two-dimensional object function has passed 
through a spatial filter to whose output random noise has been added to 

( 4 - 9 ) 

produce the observed data. On the basis of these data and the 

known filter transfer function, the original object function is to be esti- 
mated with minimum mean-square error. 

When imaging incoherent light, it is the radiance of the object 
plane that is to be estimated, not the value of the light field itself. The 
light field is a complex spatio-temporal Gaussian random process, to 
whose mutual coherence function the radiance is related. This radiance 
like its counterpart in the image plane, is a quadratic functional of the 
field. Furthermore, linear filtering theory does not explicitly distin- 
guish the two basic types of noise, that due to the inherent fluctuation of 
the light fields of object and background, and that associated with the 
process of recording the light. Instead, it makes certain assumptions 
about the statistical properties of the noise, which is usually assumed 
to be white and Gaussian and to leave such parameters as its mean- 
square amplitude and its bandwidth to be measured separately. Such a 
treatment is incapable of determining the fundamental limitations on 
restoring degraded images. 

Helstrorr/^^has shown how the statistical properties of the in- 
formation-bearing quadratic functional of the light field modify the usual 



Shannon formula^ ^ for the information transfer from an incoherently 
radiating object. We will take into account the physical nature of the 
light and the statistical properties of the recording process, and 
analyze an imaging system in which the light from the object plane is 
focused onto a photosensitive surface, from which it ejects photo - 
electrons. The surface is divided like a mosaic into a large number of 
small, insulated spots, from each of which the photocurrent can be 
measured. The measured values of the photocurrent constitute the 
data on which is to be based an estimate of the radiance distribution 
function of the object plane, which corresponds to the 'true image'. 

In chapter 1, we first describe the statistical properties of the 
light field emitted by the incoherently radiating object, and then calcu- 
late the statistics of the instantaneous illuminance in the image plane, 
which is assumed to be a measurable quantity that provides the data on 
which a linear estimate of the object radiance function will be based. 

The object radiance function is considered as the sum of a uniform 
background and one of an ensemble of spatial random processes, 

which form the class of object to be estimated. The linear filter has 
the same form as for the estimation of the signal degraded by signal- 
independent additive noise. The equivalent noise spectral density con- 
tains a contribution not only from the background white noise, but also 

from the object- light fluctuation and the uniform object background, B^. 

• • ( 5 ) 

This is very different from the model treated before. 
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In chapter 2, we will discuss the interaction of the light field 
with the recording medium. When the product WT of the bandwidth W 
of the incident light and the observation time T is much larger than 1, 
the numbers of photoelectrons emitted from the spots of the photo- 
electrically emissive recording surface are statistically independent 
with Poisson distributions whose mean values are proportional to the 
illuminance at the spots. The linear least- squares filter for the esti- 
mation of the object radiance function has the same form as those 
described by the Wiener filtering theory, but with an equivalent noise 
spectrum equal to the sum of a constant term arising from the inter- 
action of the uniform object background and the background white 

noise with the recording medium, the shot noise, and a spatial variable 

term due to the fluctuation of the light field produced by object and back- 

( 12 ) 

ground noise. In particular, it is shown that under normal circum- 

stances the shot noise far exceeds the noise associated with random 
fluctuations of the object and background noise. The shot noise is 
signal-dependent. Linear least-squares estimation, which is optimum 
for the situation when the signal is corrupted by signal- independent 
additive noise, is inadequate for the restoration of images degraded by 

( 13 ) 

signal-dependent noise. The principle of maximum likelihood leads 

to nonlinear equations for the estimates, which can be solved by linear 

(14) 

iteration or by methods of nonlinear programming. The nonlinear 

estimator is signal-dependent and is asymptotically convergent to the 
linear estimator, when the received light is strong. 
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To provide the insights about the practical implications of the 
analytic result for the maximum- likelihood estimator, computer simu- 
lated experiments are performed to demonstrate the significance of the 
nonlinear estimation scheme in chapter 3. One-dimensional objects 
distorted by motion blur and diffraction have been estimated in com- 
puter simulations. When the background illuminance is weak and the 
contrast of the object is high, the maximum- likelihood estimate is 
much superior to the linear estimate. 

In part II of this thesis we will present new methods of optical 
signal processing using only coherent optical systems. 

Optics and electronics are means of information processing 
which have many common features. Communication theory is well 
developed in electronic signal processing. Recently it has been suc- 
cessfully adapted by optical physicists, especially in the development 
of holography^ ^ and spatial filtering. Many optical phenomena 

can be easily understood from the point of view of communication 
theory. Modulation, demodulation, heterodyning, harmonics genera- 
tion, and filtering have been used successfully in electronic signal 
processing; we can find their counterparts in optics. In optics, modu- 
lation can be achieved by carrier -frequency photography^^' gen- 
eration of harmonics by grating diffraction; demodulation and hetero- 
dyning by the moire effect; and filtering by spatial filters. Carrier- 
frequency photography, diffraction by gratings, and the moire effect 
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all have their unique functions in optical signal processing. Two new 
types of interferometers have been developed from the combined use 
of those concepts. 

Diffraction by gratings and the moire effect have been used 

( 21 ) 

extensively in optics. Gratings have been used as beam splitters. 

Ronchi invented a shearing interferometer using a grating as a beam 

divider, which has been used extensively for testing optical components 

such as lenses. This technique is relatively simple and inexpensive. 

But the use of just one grating results in many drawbacks, such as 

(a) the fixed amount of shear, (b) the dependence of the fringe contrast 

on the degree of equality in light flux between different diffraction 

( 2 1 ) 

orders, and (c) the overlap of many diffraction orders. 

The drawbacks of Ronchi's interferometer can be overcome by 
the use of two superimposed gratings instead of just one grating as the 
beam divider. The amount of shear in Ronchi's interferometer is pro- 
portional to the spatial frequency of the diffraction grating used. If we 
can vary the spatial frequency, we can have variable shear. This can 
be done by the moire effect of the two superimposed gratings. When 
two gratings are superimposed, owing to the moire effect a variable 
spatial frequency is produced by rotating the angle between their axes. 
An interferometer, called the "grating interferometer, " has been 
developed from two superimposed Ronchi rulings used as beam 
splitters. All the drawbacks in Ronchi's interferometer are overcome. 
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Gar rier -fr equency photography has some similarity to holog- 
raphy. In particular, both are two-step processes, and diffraction 
methods are generally used in the retrieval stages. Holography has 

been used very successfully in interferometry. Holographic inter- 

( 22 ) 

ferometry can be operated either in real time or in two stages. 

A hologram's ability to record both the amplitude and phase of 

a wavefront is the reason for its success in interferometry. Carrier- 

frequency photography can register the amplitude as well as the phase 

variation of objects. When a grating is imaged through a complex 

object, the resultant image is an amplitude -modulated version of the 

grating with the slits shifted according to the phase variation of the 

object. This kind of encoding scheme, involving both amplitude- 

modulation and pulse -position modulation, called the "detour phase" 

( 23 , 24 ) 

by A. Lohmann, is the essence of computer -generated holography. 

It is very difficult to retrieve the information about the phase varia- 
tion of the object directly from the distorted grating image. But if 
we "beat" the distorted grating image with a regular grating (or 
master grating), moire fringes appear which will reveal the phase 
variation. A new type of interferometer, based on carrier -frequency 
photography and the moire effect, has been conceived. It is called 
the carrier -frequency interferometer. It can be operated in real 
time or in two stages. The system is achromatic and white light can 


be used. 
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In chapter 4» we briefly review diffraction by gratings, the 
moire effect, carrier-frequency photography, and achromatic holog- 
raphy. This review will serve as the foundation of the following 
chapter s. 

The grating interferometer is the subject of chapter 5, and the 
carrier-frequency interferometer is described in chapter 6. Many 
kinds of masks can be used as modulation carriers in the carrier - 
frequency interferometer. We use linear gratings, checkerboard 
gratings, and circular gratings. 

The applications of these two types of interferometers in opti- 
cal information processing are presented in both chapters 5 and 6 for 
each interferometer. New methods for optical testing, phase object 
visualization, image polarity reversal and color complementarity, 
two-dimensional photometry, spatial differentiation, complex ampli- 
tude addition and subtraction, image multiplexing, pure spatial sinu- 
soidal wave generation, fabrication of filters, focusing, and alignment 
are presented. Some of the applications can be used in optical com- 
puting, general image processing, and two-dimensional photometry. 
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Part I. Restoration of Images Degraded by Signal-Dependent 

Noise 

In this part of the thesis we will discuss the restoration of 
images degraded by noise that is neither Gaussian nor additive. 

When imaging incoherent light, it is the radiance distribution 
of the object plane that is to be estimated, not the value of the light 
field itself. The object radiance function, like the illuminance in the 
image plane, is a quadratic functional of the field. We are going to 
discuss the effect of the quadratic rectification of the field and the 
interaction of the light field with the recording medium in our imaging 
system. The optical system used is an imaging system in which the 
light from the object plane is focused onto a photosensitive surface, 
from which it ejects photoelectrons. 

We begin by describing the statistical properties of the object 
light field and the image -plane light field. A linear estimator is de- 
rived for the restoration of the object radiance function based on the 
measured data which is assumed to be the illuminance in the image 
plane. The interaction of the light field with the recording medium 
is then taken into account. Linear and maximum- likelihood esti- 
mators of the object radiance function based on the measured numbers 
of photoelectrons emitted from the recording medium have been de- 
rived. Computer simulations are performed to demonstrate the sig- 
nificance of the maximum- likelihood estimate, which is superior to the 
linear estimate when the light level is low. 


11 



Chapter 1 


Statistical Description of Object Field and Linear Estimation 
of Object Radiance Function Based on Measured 
Image Illuminance 

In this chapter, we will first discuss the statistical representa- 
tion of the object field and derive the linear least-squares estimation 
of the object radiance function based on the illuminance distribution in 
the image plane, which is assumed measurable. The interaction of 
the light field with the recording medium will be neglected here and 
will be discussed in the next chapter. 

The electromagnetic field emitted by the object is assumed to 
be linearly polarized, quasimonochromatic, and represented by its 
"analytic signal. 11 ^ The analytic signal of the object light field will 
be considered as a circular complex spatio-temporal Gaussian random 
process which is stationary, ergodic, and spectrally pure. Linear 
scalar diffraction theory is used for wave propagation. The relations 
among the covariance function of the object light field, that in the 
image plane, and the point spread function of optical system will be 
discussed. The radiance function of the incoherently radiating object 
is considered as one of an ensemble of spatial random processes, 
which form the class of the objects to be examined. The background 
noise is assumed to be a circular complex spatio-temporal, white 
Gaussian random process with spectral density N. The instantaneous 
intensity of the light field in the image plane is assumed to be 
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measurable and to constitute the data. A linear spatial filter based on 
these data is derived for the restoration of the object radiance function. 


1. 1. The Statistical Properties of the Object Light Field 

The electromagnetic field emitted by the object is assumed to 
be linearly polarized and quasimonochromatic. The light field radi- 
ated by the object at point u = U y ) of the object plane can be 
specified by a scalar function 

U n (u,t) = Rf tJu.t) e lQt , (l.D 

where Q, ~ 2ttc/\ is the central angular frequency, X the wavelength 

at the central angular frequency, and iJtq(u, t) the complex envelope 

of the field. "Rf " indicates the real part of a complex number. 

ty (u, t) e 1 ^ * is called the analytic signal. The object light field can 

be completely represented by the complex envelope of the object light 
( 2 ) 

field, t 0 {u, t). 

The complex envelope of the object light field is a circular 
complex spatio-temporal random process with the complex auto- 
covariance function 


VW W 


1 < *0 ( »l ,t l ) *0 ( S2' t 2 )> 
*0 ( W W • 


( 1 . 2 ) 


where ( • ) denotes the ensemble average, and * means the com- 
plex conjugate. The complex autocovariance function is proportional 
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to the mutual coherence function with a factor exp[iQ “ t^)] * 

Since the object is a circular complex Gaussian random process, the 

(3) 

following relations hold: 

VvV* = 0 ' (1,3) 

and 

<*0 ( Sr‘l ) ♦oWV VVV *0* ( »4' t 4 ,) 

= (♦ 0 ( Sl- t l ) * 0 *(u 2 ,t 2 )> <* 0 (S 3 .t 3 t**(V t 4 )> 

+ <* 0 ( Sj' t l ) ' t 0 ( »4' t 4 )> < +0 < »2 ,t 2 ) *0 (li 3 lt 3 )> 

= 4 WV W »o ( »3'»4 ; W 
+ 4 »0 < Hl’V t V t 4 > « 0 ( W W ' U-4) 

It can be shown that in general a quasimonochromatic wave- 

(4) 

field is spectrally pure. The spectrally pure light beam is charac- 
terized in such a way that the superposition of light beams will not 
affect the spectral distribution; the covariance function is then reduc- 
ible to the product of a spatial covariance function ) and a 

temporal covariance function \{tj - t^), i. e. , 

t l* t 2 ) = 0 O ( -l'~2 ) X (t l ' t 2 ) ’ (1 * 5) 


As a temporal random process, the light field is considered as sta- 
tionary and ergodic. The temporal covariance function is usually 
normalized so that = ^ • 
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The bandwidth W of the object light field is defined as 


w= lx(0)| 2 



dr 


■[/ 


Xfo) dew /Z 



/ 


X(iw)| daj/2n 


( 1 . 6 ) 


where 


X(«u) 


00 

f iw T . 

J X ( T ) e <1 T , 
- 00 


(1.7) 


is the Fourier transform of the temporal covariance function, or the 

temporal spectral density. The bandwidth of the object light field is 

so defined that if X((w) is constant over frequency range 

- ~ W < uu < r - W and zero otherwise, then the bandwidth defined 
4 4 

in Eq. (1.6) equals to W. ifr (u,t) and i)f^(u, t + W *) can be con- 
sidered statistically independent, and W ^ is usually called the 
correlation time. The correlation time of the signal is inversely 
proportional to its bandwidth. 

For a spatially incoherently radiating object, the mutual 
coherence function takes the form 

V V = (x2/4tt) M}Ji - u 2 ) X (tj - t 2 ) , (1.8) 

where 6 (u) is the two-dimensional delta function. B(u) is the 
radiance function of the object and can be decomposed as 
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into a known uniform background and a deviation b(u). The 

deviation b(u) represents the scene of interest and is what we 
want to estimate. The radiance is so defined that B(u)/4tt is the 
power emitted per unit area per steradian in a direction normal to the 
object plane. 

Although the object light field is not spatially stationary, the 
object radiance function is very similar to the spectral density of a 
time series. ^ To estimate b(u) is equivalent to estimating the 
spectral density in time series analysis. The variability and unbiased- 
ness always conflict with each other. ^ 


1.2. The Light Field in the Image Plane 

The system we want to study is shown in Fig. 1. 1. An optical 
observing aperture A is used to collect the light emitted by the object. 
The coordinate variables u, r and x represent the object, aperture 
and image plane respectively. The aperture contains a narrow-band 
filter used to block the background light outside the spectral region of 
the object light, and a lens to focus the light onto the image plane 1. 

The distances from the aperture to the object plane and image plane are 
Z^ and Z. respectively, and satisfy the image-forming relationship 



where f is the focal length of the lens used. 



Radiance 



0 


Figure 1. 1. The object and image -proces sing system. 

O = object plane, A = aperture, I = image 
plane. A narrowband spectral filter for 
object and background light is not shown. 
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If the object, aperture, and image planes are assumed to be 
far apart, and the aperture is not too big, then all light rays can be 
considered paraxial, and it suffices to use scalar diffraction theory to 
describe the wave propagation from the object plane to the aperture 
plane and then to the image plane. 

At point r of the aperture plane A, the total field is 

♦ a <£. t) “ * g (r. t) + * (1.10) 


where 4 (r,t) represents the background noise, and 
n ~ 

d2 ~ 


U. ID 


represents the component due to the object, or signal. Here 

S(u,r; t) is the instantaneous amplitude point spread function describ- 

A J A/ 

ing the wave propagation from the object plane to the aperture plane, 
and 0 indicates an integration over the object plane. The point 
spread function S(u, r; t) can be a time varying function, such as in 
the situation that the object light field propagated through a turbulent 
medium, or when the object is moving relatively to the aperture plane 
in the period of observation. For a diffraction limited system, in the 
paraxial approximation, it is 

S 0 (u,r; t) = S Q (u,r) 

= i(X Z 0 ) _1 exp(iKZ 0 +iK | u-r| 2 /2 Z Q ) , (1.12) 

where K = 2tt/X is the propagation constant. 
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The background noise (f (r, t) is a spatial and temporal white, 
circular -complex Gaussian random process whose mutual covariance 
function, after filtering by the temporal filter, is 


*nter£z ;t rV 2 ^n-l’V ^n^Z’V^ 


= N6( £r^ ) ^ (t l _t 2 ) ’ 


(1. 13) 


where N is the spectral density of the noise. 
The field at the image plane is 

^.(x.t) = / t a (r f t) S a (r,x) d 2 r 
A 

= i|) . (x, t) + i|i . (x, t) , 

Y si T ni ~ 


(1. 14) 


where S (r,x) is the amplitude point spread function between the 

El ^ ^ 

aperture and the image plane, and in the paraxial approximation 

S a (r,x) = i(XZ.f 1 exp[iKZ. + iK | r - x| 2 /2 Z. - iKr 2 /2F ] . (1.15) 

The last term in the exponent is the quadratic phase factor introduced 
by the lens. 

The field due to the background noise is 


'll .(x.t) 

ni 




t 


(1. 16) 


where A indicates an integration over the aperture plane. And the 
field due to the object light is 
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(x, t) = / t]r Jr, t) S (r,x) d 2 r 


Si ^ 


■ // 
A 0 


V-' t} S 0 ( -’- : S a { Z’~ ) ^ d2}i 


-/ 


i([ (u, t) S(u, x; t) d u 

T (J ^ w <v ^ 


(1. 17) 


where 


S(u, x; t) = / S (u, r; t) S (i , x) d 2 r , 

A/ #v J (J ^ a ^ ^ V 


(1. 18) 


is the amplitude point spread function between the object and image 

(13) 

planes. For a diffraction limited system 


S (u, x; t) = S(u, x) 


= U 2 z 


“»[‘f (f ’ f )] 




d 2 r . 


(1. 19) 


If the aperture is large enough, the integral in Eq. {1* 19) will 
approach a delta function located at x = -(Z./Z^)u, and the field at 
point u will form its M geometrical image n at 


x = -(Z./Z ) u 

— i o ^ 


= - Mu 


{1.20) 
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where 


M = Z./Z 
1 o 


( 1 . 21 ) 


is the magnification factor. 

It is more convenient to deal with image coordinates having the 
same scale as those of the object, and we will define a new coordinate 

V = -<Z /Z.) x 

/v o 1 ~ 

= -M~ 1 x . (1.22) 

The diffraction-limited point spread function, from Eq. (1.24), in the 
new coordinates is 


S(u, v) 


U 2 z Z f 1 e X p[iK(Z o+ Z )]e X p[i|-(u 2 + f v 2 ) 

o' o 


f exp[iK(u-v) - r/Z ] d 2 r 

1 V /y/ n J A/ 


= (X Z Z Z.r 1 exp[iK (Z +Z 
O 1 L o 


.)] exp£ 


•] 


— - (n 2 +Mv 2 ) | ♦ A • J(u-v) , 


(1.23) 


where 


J(u) = A 1 J exp(2rr i u • r / X Z ) I. (r ) d~r 
A 


(1.24) 


is the Fourier transform of the aperture indicator function 

( 1 r € A 

l.(r) = { 

( 0 r ^ A . 
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Since the noise background field ty^(r, t) an d object light field 
are independent, the covariance function of the field in the image plane, 
from Eq. (1. 14), is 


t’&l’Zz'’ W = 2 *i ( X 2 - t 2 ) > 




* 


+ 7 < i|> .(v , t ) i|i .(v , t )> . 
2 Y ni ~1 1 ni ~2 2 


(1.25) 


In the image plane the covariance function of the light field originating 
from the object, from Eq. (1. 17), is 


0 .(v.,v_; = ~ < ^.(v , t ) i])* (v , t )> 

si ^1 12 2 Si^l 1 T si~2 2 


= (X 2 /4tt) / B (u) S(u, v^ t ]L ) S (u, v 2 ; t 2 ) d 2 u • X^"^)* (1.26) 


and for the diffraction-limited case 


2 2-2 2 2 ”1 
- (X / 4tt)( X Z q ZJ X (tj " t^) exp I ~ (Vj^ ■ v 2 M 

LIZ J 


(4ttZ 2 ) 1 A exp 


A 2 f B (u) £ (u - v ) «9 (u - v ) d 2 u 


lnZ i , 2 2,1 , 

“2 - X 2 >J • x < t x - 1 2 > 


(1.27) 


j B J(v -v ) + (X Z f 2 A J b(u) J(u-vj /(u-v.) d 2 u \ . 

^ U ^ 1 ~2 O Q ** * ** ~2 ** ? 
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The background noise field A . (v, t) has covariance function 

ni 


0 .(v , v ; t , t ) = - ( f .(v ,t.) i|i (v , t )) 
m ~ 1 "2 12 2 m ~1 1 ni ~2 2 


= N( \Z.) 


-2 


exp 




in Z. 


xz. 


. 2 2 . 
(v - v ) 
^1 ~2 






(1.28) 


After combining Eqs. (1.25), (1.27), and (1.28), the covariance func- 
tion of the light field in the image plane is 




vv = 






(1.29) 


0-(v., v„ ) 

i ~1 ~2 


(4TT Z?) 


i r 

A • exp I 


r- Zi , 2 2,1 

9 ( V 1 ) 

Lz* - 1 ~ 2 J 


| B ■ ) (iti-x 2 ) + uz o r 2 A / b(« 


b (U) ^'()i-Y 2 ) d 


‘4 


where 


B = B Q + (4tt/X ) N 


(1.30) 

(1.31) 


is the equivalent background noise spectral density. Here we can see 

that the uniform object background, B q , contributes to the noise 

spectral density. This is very different from the linear model of 

signal formation, in which the noise and the object radiance function 

are independent and combined additively. In the linear model, the 

known uniform object background B q does not contribute to the noise 

(7) 


spectral density. 
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1. 3. The Instantaneous Image Plane Light Field Intensity 

In this chapter, we will assume that the instantaneous light 
field intensity is the measurable datum. A linear least-squares esti- 
mation of the object radiance function based on this measured datum 
will be derived. 

The light field in the image plane is not measurable and is not 

(8) 

of interest. The measurable datum is its intensity 

I.(v.t) = ~ H.(v,t) | 2 . (1.32) 

The mean value of the intensity is 
<I.(v,t)> = | < H.(v,t)| 2 ) 

= 0 . (v, v; t, t) 

1 -v ^ 

= (4ttZ 2 ) 1 A j B + (\ Z q ) 2 A /b(u) | J(u-v)| 2 d 2 u 
\ o 

(1.33) 

The mean value of the image plane intensity is linearly proportional to 

the object radiance function. For an incoherent imaging system, the 

illuminance in the image plane is linearly related to the object radi- 

2 

ance with point spread function K(u f v) = | S(u,v)l 

The covariance function of I . (v, t) is 

1 ** 




= 4 ^V-l’V' 2 M-2'V'^ " 4 < Ui ( Xi' t i ) l 2 > <1 

= jlO^vOvV ^ 2 

= 1 0 itei >zr Wl 2 

= (4tt zh 2 A 2 • I x (tj-t 2 )| 2 | B 2 | J(v 1 -y 2 )| 2 
+ 2BA(XZ o )‘ 2 « s (y l - v 2 ) / b(u) ^(u-y^ J"(u-v 2 ) d 2 u 

+ (XZ o )“ 4 A 2 jf b( Hl )b(u 2 ) cS^-y^ c9*{u r y 2 ) /tu^) ^-y,,) 

d2 -l d2 -2 } * (1#34) 

and its variance is 

VartMv, t) = U ^y, v; t, t) | 2 

= [ < I.(y,t)> ] 2 . (1.35) 


Based on the statistics of the measurable data, I.(v, t), a linear 
least-square estimation of the object radiance function will be derived 


in the next section. 
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1.4. Linear Least-Squares Estimation 

In this section we assume that the illuminance in the image 
plane is available for linear spatial filtering, the result is 

T 

J(u) = £ / / L(v) I . (u - v, t) d 2 v dt , (1. 36) 

10 


where I indicates an integration over the image plane, T is the 
observation interval, and L{v) is the point-spread function of the 
spatial filter. The output of the rectified and filtered image plane 
light field J (u) is 


T 

<J(u)>= // L(v) < I (u - v, t)> d 2 v dt 

10 

~f L(v)(4rr Z?) A B +(\Z ) 2 A /b(u' )|c9(u'-u-tv)t 2 d 2 u' d 2 v 

i^ 1 1 O -/ a/ I /V iv /u I ^ 

I 0 

= J B + (4rrX 2 Z 2 Z 2 )' 1 A Z j^L{v) b{uO lJ(u'-u+^| 2 dVd 2 v , 

(1.37) 


where 


J„ = (4rr Z 2 ) 1 A 
£) 1 


B 


/ L(v) d 2 v 

A/ A** 

I 


is the contribution from the background. Equation (1. 36) can be re- 
written as 

( J(u)> = J + C f A (u -u') b(u') d 2 u y , 

*+ D J Al Al 


(1. 38) 
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where 

A (u' ) = J L(v) | «9(u/ - v) j 2 d 2 v (1. 39) 

and 

c = (4 ttx 2 z 2 zV 1 a 2 . 

If we could choose a linear filter L(v) such that A(u) is a 
delta function, then the estimate^ ^ 

b(u) = C" 1 [J(u) - J B ] (1.40) 

of the deviation of the object radiance function would have a mean value 

A 

equal to b(u) itself. Then the estimate b(u) of b(u) would be un- 
biased, and the mean-square bias 

6 b = < 1 < b(u)> - b(u)| 2 > (1.41) 

would be zero. But if we go too far in that direction, the mean- 
square fluctuation 

= ( j b(u) - <b(u)> J 2 > (1.42) 

begins to increase. The mean-square bias and mean-square fluctua- 

A 

tion of the estimate b(u) conflict with each other as in temporal 
power spectral estimation, ^ We have to choose a linear filter L(v) 
which achieves the minimum of error & a + , The mean-square 

bias, calculated from Eqs. (1,37) and (1.40), is 
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e L = b 2 (u) - 2b(u) / A (u-u') b(u') d Z u # 

D ^ A / A . 


+ / / b(u ) b(u ) A (u -u' ) A(u - u' ) d 2 u' d 2 u' » 


and the mean-square fluctuation is 


£ = var (b (u)) 


= c 2 var (J (u)) 


The variance of J(u) , from Eqs. (1.36) and (1.33), is 


T T 

Var (J(u)) - ~r / J // L(v. ) L(v ) cov f I. (u-v. , t. ), I (u-v-,t )} 

~ „2 „ „ r . ~1 ~2 l ~ '*•1 1 l ~ ~2 2 


T 0 0 I I 


d 2 v d 2 v_ dt,dt 
-1 ~2 12 


2-2 2 
= (WT) (4nZ. ) A 


// L &i> L te 2 ) | b 2 1 -y 2 M 


+ 2(XZ )' 2 BAJ(v,-v,) Jb(uj «9(u-u+ Vl ) c9*(u -u +v 0 ) d 2 u, 
o ^ I ^ 1 ^ ^ 1 v j ~2 ^ l 


+ (XZ ) 4 a 2 //b(u ) b(u ) J(u -u+v ) cS*(u -u +v ) 
o ^ 1 ^ 1 — ** ^2 


X J (vu-u+v ) J(u -u+v- ) d 2 x 

^4 j[ 4 ^ -v ^ « 


(1.44) 
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Any signal processing system is designed for a particular class 

of signal. The statistical properties assigned to the class reflect the 

nature of the signal to be processed. The prior knowledge about the 

. . ( 7 , 12 ) 

statistical properties of signals is essential in signal processing. 

Here we will consider the scene of interest b(u) in the object plane as 
a realization of a homogeneous stochastic process with mean zero, 
and covariance function 

b( -l } b( ~2 )3 = "-2 ) ' (1 ' 45) 

and 

CT Z = cp (0) * (1.46) 

b b ** 


E denotes the average over the ensemble of the spatial processes 
b 

b(u). The structure of the picture, i. e. , b(u), is to some extent em- 
bodied in the covariance function cp, (u), whose width will be of the 

b ^ 

2 , 2 

order of the average size of details in the picture. The ratio /B q 
specifies the mean-square contrast. By averaging over the ensemble 
of possible pictures, we can express the mean-square error as 


e = E b { < [b(u) - b(u) ] Z > } 

= E b {<[b{u)-<b(u))f> + <[<b(u)> - b(u)f> ) 


( 1 . 47 ) 



where 

l x = E b [ < [ b(u) - <b(u)> ] 2 > = E b { £ y ) U.48) 

represents the mean-square error due to the fluctuations of the light 
field itself, and 

t 2 = E b [ ( [ k (u) “ ] 2 >3= E b {& b ) d-49) 

represents the mean-square bias arising from the wave propagation 
and the post-rectification filtering. 

From Eqs. (1.48), (1.43), and (1.45), after some calculation, 

&. = C ~ 2 • E, { Var (J(u))} 

I D ^ 

= c" 2 (WT)“ 1 (4nz 2 )“ 2 a 2 !!u Xl ) L(y 2 ) |b 2 \M Xl -z z )\ Z 

+ a 2 (U )' 4 J Jy (u -u ) c9(u -u+v ) c9*(u -u+y ) 

00 - 

x c9*(u -u +v ) J(u -u +v_) d 2 u d 2 u ? d 2 v d 2 v, 
a- 2 ^ ^ 1 ^ ^ Z 'v l ^2 | ^ 1 ~2 

= (WTr 1 (XZ o ) 4 A‘ 3 |^B 2 /|A (r)| 2 I®’(r) d Z r 

+ // 1 A( £>! 2 * b < Sl I A >( I'J ) ^ I ’ <1>50) 


where 
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$(r)=(XZ o ) 2 / cp b (u) exp(-iKu • r/Z o ) d 2 u , (1.52) 


I (2) (r) = A" 1 / I A (u) I ( r +u)d 2 u , 
A ~ •' A ~ A ~ ~ ~ 

A 


(1.53) 


I^(r, s) = A 1 / I (u) I (u +r ) I <u +s) I (u + r +s) d 2 u . 
A ~ . A ~ A'-' ~ A ~ ^ A ~ ~ ~ ~ 

A 


(1.54) 


(2) 

A(r) is the transfer function of the linear filter L{u). I (r) is the 
correlation function of the aperture function. Combining Eqs. (1. 50), 
(1. 37), (1.40), (1. 51), (1.52) and (1. 53), we find 




-I 


2 - 1 (? \ 2 2 
I 1 - <*V A ! A <£.> a <£M <1' • 


(1.55) 


and the mean square-error is 


e= e, + 


= /|A(r)| z | 


UZ q ) 4 A 2 # b ( r) [lA ){ - )]2 + (WT)"" 1 ( XZ Q ) 4 A~ 2 ^B 2 I^(r) 


+ / * b ( ^ I A >( ^^ ) d ^]| d2 


- 2 Ri 


| (X Z Q ) 2 A “ Y A(r) $b (r)I^ ) (r)d 2 r j 


+ f$, (r) d 2 r 
J b /*• -v 


(1. 56) 
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Using the variational calculus, from Eq. (1. 56), we find that the linear 
filter which gives the minimum mean-square error is 


A (r) = 


(XZ Q ) 2 

0 



f 


r € A' 


r ^ A' 


(1.57) 


where 


D(r) = 


.txiffen 2 


(WTf 1 A -1 


I 




/ 




I A )( i'5. |d2 5 


f 


(1.58) 


(2 } 

and A' is the area where 1 (r) 0 , or the so called convolved 

A ^ 

aperture. The last term in Eq. (1.58) is the noise spectrum. It is 

clear that the bigger the time bandwidth product WT and the area of 

the aperture A are, the smaller the noise is. 

In our model, we do take into account the effect of rectification, 

i. e., the measured data is the instantaneous intensity of the 

in the image plane. The linear filter A Or) that we have obtained, as 

shown in Eq. (1. 57), is different from the one using a linear model for 

(7) 

the formation of the measured data. The first term in the noise 

spectrum originated from the background noise and the uniform object 

background ; the last term is from the fluctuation of the object 

radiance function. In the linear model of image formation, only back- 

(7) 


ground noise contributes to the noise spectrum. 
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When (WT) A (cJg/B) 2 » 1 , ft(r) - (XZ q ) 2 A [^(r)] 1 for 
most of the convolved aperture A 7 , except for those parts near the 
edge of the aperture, where A(r) increases to a finite value 
(xz ) A (a /B) (WT) instead of going to infinity. Since the diffrac- 
tion- limited system, is a bandpass filter , none of the object spectrum 
outside the passband is recoverable, as can be seen from Eq* (1.57). 
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Chapter 2 


Photon-Counting Statistics and Maximum- Likelihood Estimation 
of Object Radiance Function 

In chapter 1, we have assumed that the measurable datum is 
the instantaneous intensity of the light field in the image plane. In 
reality, the instantaneous intensity of the image light field itself is not 
measurable, and it should interact with a recording medium from 
which we can take data. ^ ^ The recording medium usually used is a 
photosensitive surface, for example, a photographic film or a photo- 
multiplier tube. We will analyze an imaging system, as shown in 
Fig. 1. 1. in which the light from the object plane is focused on a photo- 
sensitive surface from which it ejects photoelectrons. The surface is 
divided like a mosaic into a large number of small, insulated spots, 
from each of which the photocurrent can be measured. These meas- 
ured values of the photocurrent during an observation interval T will 
constitute the data on which is to be based an estimate of the radiance 
distribution function of the object plane. When WT » 1, the numbers 
of photo electrons emitted from each spot of the photosensitive record- 
ing surface are statistically independent with Poisson distributions 

{2-9 

whose mean values are proportional to the illuminance at the spots. 

In this model, the noise arising from the interaction between 
the light field and the recording medium predominates in causing the 
fluctuation of the measured data. This kind of noise, which is called 
shot noise, is signal-dependent. The linear filter for the 
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least -squares estimation of the object radiance function has the same 

form as we have derived in chapter 1, Eq. (1.57), except that one more 

noise term due to the shot noise has been introduced, and in normal 

circumstance the shot noise dominates. The statistics of the measured 

data are Poisson and as the noise is signal-dependent, the linear esti- 

( 12 ) 

mator is no longer optimum. The principle of maximum likelihood 

leads to nonlinear equations for the estimates. The nonlinear esti- 
mator is signal-dependent and is asymptotically convergent to the linear 
estimator, when the received light is strong. 

2. 1. Photon- Counting Receiver Model and Its Statistics 

The photosensitive recording surface in the image plane of 
Fig. 1. 1 is shown in Fig. 2. 1. The recording surface is divided into a 
mosaic of small, insulated spots of area S. The numbers of photo- 
electrons measured at each spot constitute the data. This receiver 
may be considered as a matrix of photomultiplier tubes. 

The photon counting statistics have been previously discussed 

(3-9) 

extensively in literature. We want to present a very brief descrip- 

tion of the situation that the incident signal in the image plane i|r.(x, t) 

1 ^ 

has large time -bandwidth product, WT » 1. 

To study the statistical properties of the number of photoelectrons, 

n(Xi), emitted by each photomultiplier tube centered at x. of the image 

plane, we need to determine the probability P (n; T) that exactly n 

£i 

photoelectrons are emitted from that photomultiplier tube during an 


observation interval of T. 



Area s 



Figure 2. 1. A photo -counting receiver model. 
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It seems reasonable to assume that the probability of emission 

of a photoelectron during a given interval is statistically independent 

of the number of photoelectrons emitted previously. We assume that 

we can write the probability of the emission of a single photoelectron 

in an infinitesimal interval At , P (1 ; At) as 

x. 

~i 

P (1; At) = as I . (x. , t) At , (2. 1) 

x. 1~1 

~ l 

where 

a = n / h v , (2.2) 

and r) is the quantum efficiency of the photodetector, h is 
Planck's constant, v is the mean frequency of the incident light. The 
intensity of the light field in the image plane, P(x, t), as described in 
chapter 1, is 

I.(x, t) = — | \jf .(x, t) | 2 . (2.3) 

1 Z 1 T l * 


If we add the additional assumption of zero probability of multiple 

occurrence in the infinitesimal interval At (A t -* 0) , we can write 

(3 5 7) 

the conditional probability distribution function' 9 9 of the emission 
of n photoelectrons in an observation interval of T , 


P x <n;T|l.<x.. t » = 


~T t E (£i ; T >f exp[-E(x.; T)] 


(2.4) 


where 
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T 

E(x; T) = as f | I.(x, t)| dt 
~ 0 

T 

= j a s / | i|i . (x, t) | 2 dt , (2.5) 

L 0 1 ~ 

is the energy incident on the detector of area s centered at x, during 
an observation interval of duration T . 

The probability distribution function described in Eq. (2. 4) is 
the conditional probability distribution function, conditional on a 
particular realization I(x, t), t € [0, T] , of a random process whose 
statistics have been discussed in section 1, 3. From Eq. (2. 4), the 
mean value of the number of photoelectrons from the spot centered 
at x is 

A/ 

< n(x) ) = <E(x;T)> 


T 



= asT 0.(x, x) , (2.6) 

1 A/ 

where 0.(x,,x ) is spatial covariance function of the image plane light 
field, given in Eq. (1.30). 

The covariance function of the number of photoelectrons 

emitted from two different spots can be evaluated as in the semi- 

(13) 

classical analysis of the Hanbury Brown-Twiss effect. We use the 
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equality of the variance of the Poisson distribution of n(x^) with its 

mean value. The correlation function of n(x.) and n(x ) is 

~ 1 “-J 

< n(x.) n(x.)) 

= < E(x.; T) E(x.; T)> + ( E(x.; T)> 6. . , (2-7) 

where 5.. is the Kronecker delta. Combined with Eqs. (2.5), (2.6) 
and (2.7), the resulting covariance function is 

Cov f n(x. ), n(x. ) } 

= ( n(x. ) n(x )> - < n(x .)> < n(x ) > 

/Vj "'■'J 

= cov { E(x.; T), E(x.,' T) j + ( E(x. ; T) } 6.. 


T T 

= (as) 2 / / |0.{x. x ; t 1 ,t 2 )l 2 dt 1 dt 2 +aaT0 i (x.,x i ) 6 
0 0 J J 


T T 


= (as) 2 | 0 .(x »x )| 2 / / 1 X <t - 1 ) | 2 at 
1 1 J oo 


dt +asT 0.(x.,x.) 5.. 
2 i ij 


(as) 2 TW 1 j0.(x., x )J Z +as T 0 (x., x.) 6.. , 

7 1 '''l J ’ 1 ~1 1 1} 1 


(2 -. 8 ) 


where W is the bandwidth of the object light field as defined in Eq. (1.6). 
Comparing Eqs. (2. 8) and (1. 34), we find that the last term in 
Eq. (2. 8) is an extra term arising from the interaction of the light 
field with the recording medium. The ratio of the first term to the 
second term in Eq. (2. 8) is 
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(as) W " 1 ( 0 i (x. > x.)] 1 

£ ^7 0 .(x.,x.) 

(2.9) 

« a s W _1 (X 2 /4rr) B 
= B T] s /4tt Whv « 1 . 

- ( 1 ) 

If we consider a scene illuminated by moonlight on a clear night, 

-3 - 2-1 o 

we have B/W equal to 1.6x10 Wm jj atX = 5150A. The 

- 13 

ratio has the value 0,8 X 10 T)s which will be further diminished by 

the reflectance of the scene* For illumination by full sunlight at the 

same wavelength, the ratio is larger than that of moonlight by factor 

6 

of 10 . The shot noise arising from the interaction of the light field, 

with the recording medium is thus under most circumstances pre- 
dominant in the covariance of n(x). 

To find the unconditional probability distribution function of 

n(x.), we have to calculate the ensemble average of P (n; T Jl. (x. f t)) 

1 ~i 

with respect to E(x,; T)* The mean value and variance of E(x; T), 

*vi 

can be calculated from Eqs. (1.36) and (1.38), are 

<E(x;T)> = asT <I.(x, t)> , (2.10) 

and 

Var [E(x; T)] = (a s ) 2 TW * 1 [< I. (x, t)> ] 2 . (2.11) 

The ratio of the variance to the square of the mean value of E(x; T) is 
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Var [E(x ; T) ] . 

— = (WT) 1 . (2.12) 

[ (E(x ; T) > ] 

If the time-bandwidth product WT is very large and ( E(x; T) ) is 
finite, then we can say that the random variable E(x; T) is convergent 
to its mean in distribution, and 

P(E(x; T)) -» 6 (E (x; T) - < E(x; T)> ) . (2.13) 

(14 15) 

This can be explained by the law of large numbers 9 ; E(x ; T) 

can be considered as a sum of WT independent identical random 
variables. (16, 1?) 

From the above argument, we conclude that the unconditional 
probability distribution is , by (2.4), 

P (n; T) = P (n; T| < I. (x; t)> ) 

X X ’ 1 ~ 

= ~ [< E(x; T)> ] n exp[-<E(x; T)> ] , (2.14) 

if WT » 1 . 

Under the general circumstances the numbers of photoelectrons 
emitted from two different spots of the photosensitive recording medi- 
um are statistically independent and Poisson distributed, if we can 
neglect the classical fluctuation noise. It is clear from Eqs. (2. 8) and 
(2. 9) that the classical fluctuation noise is negligible. So the condi- 
tional joint probability distribution function of the numbers of photo- 
electrons n^, n^, . . . , n^ emitted from the spots centered at 
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L . f X„ » » • * » 


n 2' • • * n M’ ^ I \ { -V \^>2 * * ' * ' V*M' ^ 


M 


- n 


P (n. ; T 1 I.(x. ; t)) . 

. , x. 1 1 i ~i 

J = 1 ~i 


(2. 15) 


For WT » 1 , E(x.; T) converges to ( E(x. ; T)) in probability, 

1 ~1 

and the joint distribution of E{x. ; T) , i = 1, . , M, will converge to 

~i 

M 

Tf (E(x.; T)) in probability, so that 
i = l ~ 1 


M 


P(n lf n 2 , . 


. .,n..;T)=Tf P (n. ; Tj < I.(x. ; t)> ) 
M . x. j 1 i ^-i 

J = 1 ~J 


M 


-fl 


— [<E(x.;T)>] n exp[- < E(x ; T)>] . (2. 16) 
„ „ n. ! 

3 = i i 


The numbers of photoelectrons n. emitted from the spots x^ are 
statistically independent and Poisson distributed and have mean values 


<»i> = 


<E(x. ;T)> , 

^1 


(2. 17) 


if WT » 1 . 
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2,2. Linear Least-Squares Estimation of Object Radiance Function 
Based on the Measured Numbers of Photoelectrons 

We have derived a linear estimator for the object radiance func 

tion based on the measured image plane intensity in chapter 1. In this 

section we are going to derive a linear least-squares estimator based 

on the statistical properties of the measured data, the numbers of 

photo electrons emitted from the spots of the image-plane recording 

surface. ^ 

Just as we have divided the photosensitive recording surface 
into a mosaic of spots, we sample the object plane at a set of points 
u = (u , u. ), and the radiance b(u.) is to be estimated as a linear 

-3 j x jy 

combination of the numbers of photoelectrons received. 


b. 

J 


b(u.) 


2 


L. n . + C 
jm m o 


(2. 18) 


where C is a known constant, and the weighting factors L intro- 
o J m 

duced by the processing matrix L are selected to minimize the mean- 
square error 

£ = < [ b(u) - b(u) ] Z > 

= S E b < [kfaj) - b (u.)f> * (2.19) 

j 

where E, denotes the ensemble average over random process b(u) 
b ~ 

The mean value of the estimator, from Eq. (2. 18), {2. 6), and 


(1.30), is 
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(b (u. ) > = y L. <n > + C 
~i im m 

i m J 


C + asT 0.(u.,u ) 

O 1 '-J 


C +as T 
o 


(4ttZ^)' 1 a|b + (\Z o )' Z A / b(v) 1 c9(v-U m )| 2 d Z v J 


= Cs L, fb(v) 1*9 (v -u ) 1 2 d 2 v , 
Z-t jm-' ~ 1 ~ ~m 1 ~ 

m 


( 2 . 20 ) 


where 


C = aT (4rr Z 2 )" 1 (\ Z q )“ 2 A 2 , 


( 2 . 21 ) 


and 


C = -a s T A (4TT Z 


f)-‘B I 


m 




( 2 . 22 ) 


We assume that the area of each element in the image plane s is so 
small, and the spots are so close together, that the summation in 
Eqs. (2.20), and (2.22) can be approximated by an integral. 


(b(u) ) = C JJ L>(u - w) | cS (w, v)| 2 b(v) d 2 v d 2 w , 

and 

C = -aT A(4 tt Z 2 )" 1 B / L(v) d 2 v 

O 1 J 

= -aTJ 6 . 


Equation (2.21) is identical to Eq. (1.40), except for a constant factor. 
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The variance of b(u.) can be evaluated from Eqs. (2. 18), (2. 8), 


and (1.30), and is 


Var b(u.) = / / L L. cov { n , n } 
~i 4^ if im i m 

J l m J J 


= J Yl. L. (as) 2 |0.(x, ,x )| Z T/W + /L 2 asT 0 (x x ) 
£-> ii im 1 i -'-m 1 “ if i 

Zm J " 


2 2 2 

a (T/W)//L(u.-v 1 )L(u.-V 2 )|0.(v I ,v 2 )| dvjdv 2 


r 2 2 

+ aT / L (u. - v) 0 . (v, v) d v 


(2.23) 


in the limit s -+ 0 . Comparing Eq. (2.23) and (1. 49), we know that 
the last term in Eq. (2.23) is the extra term in the variance of the 
estimator when the measured data are the numbers of photoelectrons 
emitted by the spots of the photosensitive recording surface through 
the interaction of the light incident on the image plane with the record- 
ing medium. 

Applying the same analysis as we have used in chapter 1, we 
find the linear filter which minimizes the mean-square error is 

A(r ) = C~ 1 (X Z ) 2 A (r) I (2) (r) D^fr) , (2.24) 

^ o b ~ A ~ 1 ~ 

where 

D.(r) = D(r) +4 ttB (X 2 aT A) _1 . (2.25) 

Jl ^ ^ 


and 
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Dfe) =» b te)[I^ > (£)] 2 +A' 1 (WTr 1 [B 2 I^ ) (r)+/ 4b (s)I^ 4) (r,8)d 2 sj . 

The functions I^(r) and I^(r, s) are defined by Eqs. (1.53) and (1. 54) 

A ~ A ~ ~ 

and A is the area of the aperture. The optical system that we have 
discussed here is assumed to be a diffraction-limited system. 

Comparing Eqs. (2.24) and (1. 57), we find that the only differ- 
ence between those two linear spatial filters is the shot noise arising 
from the stochastic nature of the photoelectric emission. Under most 
circumstances, as we have discussed in section 2. 1, the shot noise 
arising from the quantum nature of light and its interaction with 
recording medium is predominant in the noise spectral density of the 
linear filter in Eq. (2,24), so these two linear filters in Eq. (2.24) 
and (1. 57) are basically very different. Equation (2.24) can be 
approximated as 

Mr) =C' 1 (XZ o ) 2 A* b (r) I®\r) {® b <r )[ (r )1 2 +4 ttB( \ 2 a TA)‘ 1 ] 

(2 \ - 1 — 2 - 1 - 1 f? \ _? 

= C 2 [I A ( ~ } ] £ 1+4ttB ^ aTA > $ b (r)[I^ ; (r)] }. (2.26) 

The signal-to -noise ratio for this system is 
Y = X 2 aT(4 1 TB)" 1 $ b (r)[I^ ) (r)] 2 . 


(2.27) 
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2. 3* Maximum- Likelihood Estimation 

Under most circumstances the shot noise predominates in 
causing the fluctuations of the numbers m of photoelectrons emitted 
from the spots centered at of the recording medium; the noise 

arising from the fluctuation of the light field can be neglected. The 
object that we discussed radiates incoherently; it is the radiance dis- 
tribution of the object that is of interest, not the light field it radiates. 

If the noise arising from the fluctuation of the light field can be neg- 
lected, then we can deal directly with the object radiance instead of the 
field which the object radiates. The transmission of the object radi- 
ance function to the image plane can be described as incoherent image 
formation. The illuminance in the image plane is linearly related to 
the radiance distribution in the object plane by 

I(x) = f K(x, y)B(y)d Z y , (2.28) 

O 

where K(x) is the incoherent point-spread function. For a diffraction- 
limited optical system 

K(x-y) = 1 S(x, y) | Z , (2.29) 

A* ->u< 'W ^ . 

where S(x, y) is the coherent point-spread function and is defined in 
Eq. (1.23). 

For large time -bandwidth product, when the fluctuation of the 
light field can be neglected, the numbers are statistically inde- 
pendent, discrete random variables with Poisson distributions and 
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mean values given by 


< n. 1 B(y)> = asT I(x . ) 
1 ^ 1 


(2. 30) 


as has been shown in section 2. 1. Under these Poisson statistics, the 
shot noise predominates in causing random variation of the observed 
data. The shot noise is signal-dependent. Linear least-squares esti- 
mation, which is optimum for the situation when the signal is corrupted 
by signal-independent additive Gaussian noise, is inadequate for the 
restoration of images degraded by non-Gaussian-distributed, signal- 

dependent noise. ^ We will derive a nonlinear estimation scheme 

f 1 8 19) 

based on the principle of maximum likelihood. ’ A nonlinear 
estimator has to be solved iteratively. The estimator is signal- 
dependent and varies in each step of the iteration. A digital computer 
will be used, and we will derive our estimator in digital fashion. 

Sampling the object plane at M points U ^ , we replace the 
integral by a sum and write for the mean numbers of photoelectrons, 


<n. 1 (b.}> = asT (I o + V ) 


where 


V. 

i 


m 


I 


K..b. , i 

U J 


1 , 2 , 


N 


b. 

J 


b(u.) 

M 

I 

j = l 


K.. 

U 


( 2 . 31 ) 


( 2 . 32 ) 


I = B 

o o 


v i . 


( 2 . 33 ) 
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Here I is a uniform, known illuminance, to which background light 
may contribute, the b/s are the quantities to be estimated on the 
basis of the observed numbers n. of photoelectrons, and K_ is a 
point-spread matrix derived from the point-spread function K(x, y). 
Equations (2. 31) and (2. 32) can be written in matrix form 

<n|b> = asT (I o 1 +V ) , (2.34) 

and 

V = K b , (2. 35) 

where n is a (N x 1 ) column vector with elements n. , V is a 

~ l ~ 

(Nx 1) column vector with elements V\ , K is a (Nx M) matrix with 

elements K.. and 1 is the column vector of l's. 

!J 

The joint probability distribution function of n.'s is 

n. 

N [asT (I +V.)] 1 

P^nJIfb.}) =J\ — exp[- a sT<I o +V.)]. {2.36) 

i = l r 

This joint probability distribution function is conditional on a particu- 
lar realization b(u), which is a spatially stationary random process 

with zero mean value and auto covariance function <b (u) . The un- 

b 

conditional mean value vector and covariance matrix are 
<n> = E < n j b > = asT I 1 , 


and 


(2.37) 
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<fc- <i}>> <5- <n» T > 

= (asT)^ K 0 K T -f asT I I , (2.38) 

a, ^ O 

where I is the identity matrix, and 0 is the covariance matrix 

0 b = H 0 ij I with 

0 ij = 0 b ( -i ‘ V ’ (2 * 39) 

the covariance of the radiance at point u. and u. . From Eq. (2. 38), 

^-1 

we can see that the first term in the right-hand side is related to the 
signal and the second term is related to the noise. The signal-to-noise 
ratio is equivalent to 

R c = asT (a* /I ) . (2.40) 

o bo 

A 

The linear estimate, b , of b based on the measured data n 
can be easily derived from Eqs. (2. 37) and (2. 38) and is 


b 



T 

K K + 


I 

O 

asT 




K 


(J- 

\ asT 


n 



(2.41) 


The linear estimate b given by Eq. (2.41) is identical to the one given 
by Eq. (2.24) if the noise arising from the light field fluctuation can be 
neglected. 

In using the method of maximum likelihood, which requires the 
knowledge of the probability distribution of the estimate, we have 
already assumed in applying linear estimation that the object radiance 
function B(u) is a spatially stationary random process with mean value 
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B and autocovariance function 0 , (u) . Here for simplicity we will 
o b ^ 

assume that it is Gaussian distributed. The width of 0 ^(u) as a func- 
tion of u represents the size of typical finer details in the object plane, 

and the ratio 0, (0)/B specifies a mean-square contrast. The size 
b o 

of the typical finer details and the average contrast are the most funda- 
mental (and necessary) information to specify a scene. We are not 
asserting that each object radiance distribution to be estimated will 
actually look like the chaotic function by which we commonly picture a 
Gaussian process. Our assumption is rather that a method that works 
well for such an ensemble of objects will work well for most objects of 
similar structure, and structure can seldom be specified in terms 

less primitive than the size of typical details and the mean-square 

d.20,21) 


contrast. 


( 22 ) 


The joint probability density function of the estimanda b. 


is then 


M M 


P((b.)) - M o exp(-i I i Uyb.b. ], 


(2. 42) 


where M is a normalization constant and the m .. are elements of 
o M ij 


a matrix 


M ” 0 


-1 


(2.43) 


inverse to the covariance matrix 0 . 

b 

In the absence of a scene in the object plane, b(u) = 0, the 
illuminance in the image plane is the uniform illuminance, and the joint 
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probability distribution function of the numbers n. of photoelectrons is 

n, 

N (otsT I ) 1 

P (fn.}) = Tf r- 5 exp(-asTI ) . (2.44) 

o ' i J l n. ! o 

i = i 1 


The likelihood ratio of the numbers n. , i. e. , the quotient of 
their joint probabilities in the presence and absence of an image, is 
then 


A (fm) | {bj ) 



(2.45) 


The maximum-likelihood estimate of the true image is given by the 

set of values b. for which 
i 


A({n.}) = A({n.} \ {b.}) P((b.3) 


(2.46) 


is maximum. Maximizing A({n^}) is the same as maximizing its 
logarithm, since the logarithm is a monotonic function. The maximum- 
likelihood estimates of the b.’s are the set of values which minimize 


F({b.]) = -In [A({n.3)] 


N 


M M 


= - 1 [ n j *»a„ +VJ -asx a 0 +v. ) ] + | y X u jk •>, \ • 
i-l j=lk=l JJ 


(2.47) 


from which terms independent of the estimanda { b. } have been 

j 

dropped. 
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The set of values [b. } which minimize the negative logarithm 

j 

of the likelihood ratio, F( {b.} ) , can be obtained by solving the set of 

j 

equations resulting from setting the derivatives of F with respect to 

each b. to zero. Differentiating F with respect to b. , we have 
3 3 


(VF). 

J 


9F 

9b. 

3 



(2.48) 


Solving these equations for b^ 


we obtain 


b 

m 


M N 

I S* .k.. 

j=i i=i i j 


n. 


+ v. 

1 


- as T 



where 


V. 

i 


2 


3 = 1 


k.. b. . 
13 J 


Equation (2,49) can be rewritten by using 


1 


I + V. 

o 1 


I 

o 





+ v.) 

l 


(2.49) 


and becomes 
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M 

2 

3 = 1 


N 


n. 


5. + 2 k'.k.. 

jm ." 1 mi ij j2 


b. 

3 


- 2 

i = l ' i=l ml I 2 (I + 


(I +V.) 

O O 1 


(2. 50) 


where 


M 


k ' . 
mi 


- 2 


0 • k. * • 

j~l mJ 1J 


(2.51) 


We can write Eq. (2. 50) in matrix form, 


a -1 - 1 

b = A Y + A B(b) „ 

A' /\_v /V 


(2.52) 


where 


A is a M x M matrix with elements 


a 

mj 


= 6. + 


N 

2 

i = 1 


n. 


k ' . k. . ~ 
mi ij 2 


(2. 53) 


Y is a Mx 1 column vector with elements 


- 2 k'. , 

m mi \ I / 

i=l ' o 1 


(2. 54) 


and 


B(b) is a Mx 1 column vector with elements 


N 


B 


m 


■ 2 


i = 1 


k # . 
mi 


n. V. 
l i 


I (I +V.) 
o o i 


(2.55) 
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Examining matrices A, B, and Y , we find both A and Y matrices 
whose elements depend on the data but not the estimanda, and the B 
matrix depends on the received data and the estimanda. Matrix Y 
can be written as 

Y = K T / — n - ctsT 1 ] , (2.55) 

~ £ b ^ ~ ~ J 

and is very similar to part of the linear estimator in Eq. (2.41). 
Equation (2. 52) can be solved by a perturbation method. The pro- 
cedure starts with the linear least-squares estimate of the b.'s, 
which is given by Eq, (2. 41). Matrices A and Y are fixed for 
given the measured data; only matrix B is changing in each step of 
iteration. 


2.4. Comparison of Maximum-Likelihood Estimation and Linear 
Least-Squares Estimation 

The linear least-squares estimator given by Eq. (1.41) is an 
estimator independent of the estimanda, while the nonlinear estimator 
is a function of the observed data and the estimanda because of the 
Poisson distribution of the observed data. 

Equation (2.49) can be rewritten by using 


1 


I + V. 
o 1 



V. 

l 


I (I 

o o 


+ V.) 

i 


and becomes 
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b + 
m 


N 

2 
i = 1 


n. V. 

c' - 1 1 

mi I (I +V.) 
o o 1 


*5 k mi(r-° sT ) • 

i=l ' o ' 


(2.56) 


If I + V. is very large, then by the law of large numbers, 

n./U^ +V ) will converge to as T , and Eq. (2.56) can be approximated 

as 


N 

b + y k 7 . 

m *—> . mi 
i = l 


asTV i y /y \ 

i — = .2 k ^(r- asT • 

o 1=1 ' o ’ 


(2.57) 


and can be reduced to 


T 

K K 


o 

as T 


= 0 K T ( n - I 1 ) 
/ ~ ~b ^ \ CsT ~ o ~ / 


which is the linear least- square estimate. So for a strong object 
radiance, the nonlinear estimator will converge to the linear estimator. 
This can be explained by the fact that the Poisson distribution with large 
mean value can be well approximated by a Gaussian distribution. 


2. 5. The Method of Steepest Descents 

To find the minimum of F({b i ] ) as a function of {b. } , the 
method of steepest descents can be used. We will treat the estimanda 
{b. } as coordinates of a point in an M-dimensional space; points for 

J 

which F({b. }) is constant make up a level surface of the function to 


be minimized. 
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Starting at a trial point {b. } , one will get to the bottom by 

J 

always moving in the direction of the steepest slope. This direction is 
specified by the gradient 


(VF). 

J 


8F 

8b. 

J 


(2. 58) 


For point { b. } to be a minimum (or maximum) of F({b. }) , (VF). 

J J 3 

has to be equal to zero for all components, as is shown in Eq. (2.48). 

The component (VF). is equal to the 'error' by which the j-th equation 

in Eq. (2.48) remains unsatisfied. 

In the neighborhood of the trial point, F({b.}) is expanded as a 

J 

multivariate Taylor series, carried through quadratic terms. The 
next trial point is then 


b 1 = b° - X (VF) , 

J J J 


j 1*2,..., M , 


(2. 59) 


where the number X is chosen to minimize the quadratic form 


M 


M M 


([b *} ) = F((b°} ) + ^ e i (b i“ b i )+ 2 2 S C ii (b i -b f ^h^bf) 
J J i = 1 1 1 1 L i = 1 j = 1 1J 1 1 J 3 


(2. 60) 


approximating F in the neighborhood of b f where 


e, 

i 



(VF). 


and 
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c. . 
ij 


3 2 F 

3b. 3b. 

1 Jl 


b =b° 

rJ /v 


The X minimizing Eq. (2.60) is 


M 


■ I 


i = 1 


M M 

I 2 

i=l j=l 


e. c. . e. 
i ij J 


( 2 . 61 ) 


This procedure is repeated until the values of b cease to change 
significantly. 

The linear least-squares estimate is usually used as the initial 

starting point. As radiance is a positive quantity, we want 

B q + b. 0 , Vi, and at each stage of the iteration we therefore 

replace any estimate b. that is less than -B by -B + e , where 

l o 3 o 

6 is a very small positive constant. 
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Chapter 3 


Restoration of Images Based on the Measured Photoelectron 
Counts - Computer Simulation 

The maximum- likelihood estimator of the object radiance func- 
tion based on the measured numbers n of photoelectrons is highly 
nonlinear. Because of the lack of appropriate mathematical tools, it 
is difficult to analyze and thus appreciate the nonlinear estimation 
scheme directly from the estimator itself as shown in Eqs. (2.49) and 
(2.52). As a result of this that solution tends to be highly specialized. 
Both the nature of the signal and the nature of the noise affect the 
validity and usefulness of the estimator. We have tried to make our 
nonlinear estimator valid for a broad class of signals by specifying 
only its average contrast, the size of typical finer details, and its 
Gaussian probability distribution. The significance of the nonlinear 
estimator could be best demonstrated by computer simulations, which 
were performed on the Burroughs 6700 computer. Two particular 
optical system functions considered were spatially uniform motion blur 
and diffraction- limited image formation. The effects of various parame- 
ters were considered. 

3. 1. Generation of the Object Radiance Function and the Image 
Plane Illuminance Distribution 

In the last chapter, in the derivation of both the linear and 
maximum-likelihood estimators, we have assumed that the object 


61 



62 


radiance function b(u) is a stationary spatial Gaussian random process 

with mean zero and covariance function 0, (u). For mathematical sim- 

b ~ 

plicity, we will assume that the object radiance function is generated by 

( 1 ) 


a first-order autoregressive process, 


b i + i = rb i +a i+i • 


i — 1,2,..., M — 1 , 


( 3 . 1 ) 


where { z. } are mutually independent Gaussian random variables with 

2 

zero mean and variance O . The object has been assumed to be one- 
dimensional. The constant r is the correlation coefficient of two 
neighboring b.'s and has absolute value less than 1. It determines the 
size of typical details of the object radiance function. The larger r 
is, the smoother the object radiance function. The resulting process 
f b. ] has a covariance matrix 0, with elements 

i v b 


0 


.. = o_ 2 , 




( 3 . 2 ) 


where 


( 3 . 3 ) 


is the variance of the radiance function* cr^ represents the power of 
the signal. 

To generate identically distributed independent Gaussian random 
variables { , we first generate two sets of independent random vari- 
ables { } and [ uniformly distributed between zero and one. 

By using the random variable transformation 
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z. = O . f-2 In (R,.)] 1/2 cos(2ttR,. ) , (3.4) 

1 L lx £1 

the set of values f z. } will be identically distributed independent 

• l 

2 (2) 

Gaussian random variables with zero mean and variance o . 

After we generate the set of values { b. ] , we add the uniform 

background B q to it to form the total radiance function { } = 

fB +b 1 . The positivity constraint on { B . } applies here; negative 
1 o i 1 

values of B. will be replaced by zero. In other words, if b. < "® 0 > 

then we set b. = -B . A typical object { B. } is shown as the heavy 
1 o x 

line in Fig. 3. 1. Since this hard-limited nonlinear transformation has 
been applied in the generation of the scene { b^ } , the object is not 
necessarily a stationary Gaussian random process as described by 
Eq. (3. 1). However, our assumption that the object radiance function 
is a stationary Gaussian random process is mainly for mathematical 
tractability. It is most accurate when the contrast is low. The im- 
portant parameters related to the object radiance function are r, , 

and B ; r defines the size of typical details of the object, B defines 
o o 

the uniform background radiance, and defines the mean deviation 

2 2 ' 2 

of the radiance from its mean value B . The value of a, 6 / (B + cr, ) 

o bob 

is the mean-square contrast of the object radiance function, and 
2 

(asT)O b /B is the signal-to-noise ratio as shown in Eq. (2.40). 

After we have generated the true image vector B, we can 
calculate the degraded image 



True imoge 


Distorted imoge 


Original object (thick line) and its noisy version (thii 
line) as recorded after blurring by relative motion. 
The object and data are generated by using random 
number initiators (347, 5412, 3298, 9547, 2486) and 
parameters (as T) = 10, I = 3, a = 25 and r = 0. 6. 
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I. = I 1 +V 

^ 1 O ~ 

= I 1 + K b (3. 5) 

o ~ ^ 


Two incoherent point-spread functions will be used, one for motion 
blur and the other for a diffraction-limited imaging system. 

The motion-blur image degradation arises from the relative 
motion of the object plane and the image recording plane during the 
recording process. If the imaging system is perfect, and the relative 
distance of motion between the image plane and the object plane during 
the interval of recording time T is A , then the incoherent point- 
spread function is^ 


K (u) = 



0 £ u £ A 
otherwise 


(3.6) 


Here we have assumed that the relative motion is of constant speed. 
In the digital system, if A is equivalent to D sampling points, then 
the incoherent point-spread matrix can be rewritten as 


1 

1 


K 



X 


1 


k 



o 



M 



1 

1 


1J 

— >1 


T 

D 

1 


N = M +D - 1 


(3.7) 
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M is the number of points in the object plane and N is the number of 
points in the image plane. 

We have assumed that the image recording plane is large enough, 
and N = M + D - 1, so that all the light emitted by the object has been 
available for restoration. A close examination of the matrix K has 
shown that the sums / K .. for different indices i are not all equal 

T 1J 

to each other. This seems to contradict the assumption shown in 
Eq. (2. 33). The fact of the matter is that we really have to extend the 


K matrix into 


K n x 

+* D 


1 1 ... 1 
1 1 ’ . . . 1 


X 



\ 

\ 

X 


\ 


\ 






N 


x 


\ 

\ N 


X 


\ 


1 1 * ♦ • 

M + 2 D - 2 


N = M +D - 1 


J * 




(3.8) 


and the vector b has to be modified into a (M +2D -2) dimensional 
vector b' with elements 


( b . 


b/ = 

i 


i +D- 1 


l 0 


D £ i £ M + D - 1 
otherwise 


(3.9) 


The matrix K' has ^ K.. equal to a constant for all different indices 

j 1J 

i« Since we have the prior information that b/ = 0 for 1 s i s D-l f 
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and M + Dfi i i M + 2D - 2, the whole system can be fully represented 
by the vector b and the matrix K , 

/Si /V 

When the optical imaging system is diffraction-limited, the 

(4» 5 ) 

incoherent point-spread function ’ is 

(sin [W (x-y)] i 2 

K(x-y) = c | > , (3.10) 

where c is a constant, and W q is the effective width of the aperture 
of the optical system, or the spatial bandwidth of the optical system. 

In matrix form K(x - y) can be written as 


K . . 
ij 


K .. 

Ji 



= c 



[W Q (i-j)1 
i - j 


(3.11) 


In the computer simulation we have truncated K. .at some integer 
so that K. = 0 for | j| > J ; then the matrix has the form 


J, 
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K 


e 0 e l 


e l e 0 e i 

\ \ N 


\ 


\ 


\ 


\ 


\ 


\ 


\ 


\ 


\ \ 


\ 


\ 


\ 


\ 


W\ 


\ 


\ \ 

. \ * 


\ 




(3. 12) 


(3. 13) 


(3. 14) 


is the normalization constant. 

The constant is the spatial bandwidth of the diffraction- 

limited optical system. The value of 4 /W q is equivalent to the number 
of cells over which the light from each object point will be spread in the 
image plane. For example, if W - 0. 8, then 4/W^ = 5. The elements 

K. . of matrix K are (e rt , e, = e . , e_ = e _ , . . . ) = (0.268, 0..215, 

ij ** 0 1 -1 Z -4 

0. 104, 0. 021, ... ) . We can see that only e g * * e i’ e 2 ’ e 2 

are of significance. For W = 1. 2, 4 /W = 3. 33, and (e^, e^, e^, . ..) 

= (0, 395, 0.238, 0,031, . . . ), so that only e^, e^, e ^ are of signifi- 


cance. 
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A typical input object and its image distorted by a diffraction- 
limited optical system are shown in Fig. 3.1. Note that the small 
sharp peak between two big peaks has been completely smeared out. 

If we assumed that the measurable data are the instantaneous 
illuminances at points in the image plane, the linear filter for the 
estimation of the object radiance function is shown in Eq. (1. 57). It 
has been shown that none of the object spectrum outside the passband 
is recoverable from Eq. (1. 57), and the small peak in Fig. 3. 1 is not 
restorable, as will be shown in the computer experiment results. 

The argument for modifiying the vector b and the point-spread 
matrix K in the case of motion blur in order to keep I K.. equal to 

- _ .i, 

a constant for all different indices i can be extended to the situation 
of a diffraction-limited optical system. In motion-blur image degrad- 
ation, the light from the object plane has been all used in the restora- 
tion of the object, but in diffraction- limited image degradation, some 
light emitted by the object has been wasted. From Eq. (3. 12) we can 
see that the light emitted by b^ has been wasted, but all the light 

emitted by b has been collected. 

J 

3,2. Poisson-Distributed Random Number Generation 

The numbers n. of photoelectrons emitted from the spots of 
the photoelectrically emissive surface are statistically independent and 
Poisson distributed with mean values X . = asT (I + V. ) . 
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After we have generated the image illuminance function 
{ I + V\ } , we specify the parameter as T and then start to generate 
the Poisson-distributed integers m having mean values X .. . The 
probability of receiving m photoelectrons is 


X m -X. 


Pr [m X •] = — “ e 1 
1 i J m ! 


(3. 15) 


Let { Y } be a set of independent random variables uniformly distri- 
buted over the interval [0, 1] . The number of photoelectrons emitted 
by the photomultiplier tube centered at point is n. whenever 


P £ Y. < P . . , 

n. i n. + 1 

i i 


(3. 16) 


where 


n - 1 


P n = Pr{ j < n| X] = e' 


I 


(3. 17) 


j=0 


is the cumulative Poisson distribution with mean X . The probability 
that will lie in the interval described by Eq. (3. 16), and hence n. 
photoelectrons will be emitted, is precisely given by Eq. {3. 15). 

For computation purposes, we rewrite Eq. (3. 16) as 


n.- 1 
l 


Q „.-l = ^ Y i e 1 < I X i ! = Q n • (3.18) 

1 1-0 j =0 i 

x i 

We generate the set of numbers f Y. ] and then form { Y. e } , and 


n. 

1 
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begin to generate the sum { Q } , starting with n = 0, Q q = 1, and 

i 

adding each time \ n /nl . The increment X /n! can be generated 

recursively by multiplying the previous value by a factor X/(n +1). 

As soon as Q exceeds Y. e , the computer registers the number n 
n i 

as the number of photoelectrons received. 

When the mean value X . of the Poisson distribution is large, 
the Poisson distribution can be well approximated by a discrete ver- 
sion of a Gaussian distribution with mean value X . and variance X . . 
Therefore, n. can be generated by generating a random variable z. 
with zero mean and variance X ^ , and then calculating X ^ + and 
rounding it off to the nearest integer. The Gaussian random variable 
z. can be generated by the method described in Eq. (3.4). Two sets 
of independent random variables [R*] and [ uniformly dis- 

tributed between zero and one were generated and substituted into 
Eq. (3.4) to generate . We have used X = 25 as the threshold for 
the use of Gaussian distribution. 

The generated data, reduced to the scale of the image by 
dividing by (a sT) are shown as the thin line in Fig. 3.2. 

3. 3. Methods of Nonlinear Programming 

The equations for the maximum- likelihood estimate of the object 
radiance function { b. } are highly nonlinear with linear constraints. 

We can either solve Eq. (2.49) (equivalently Eq. (2.52)) or minimize 
Eq. (2. 47) to obtain the estimates { b^ } of { b^ } . The set of values 
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Figure 3.2. Original object (solid line) and its degraded image (dashed line) after 
distortion by a diffraction-limited optical system. The object and 
data are generated by using random number initiators (5431, 3214, 
5438, 3654, 2495) and other parameters (as T) = 10, I = 3, a = 15 
and r = 0. 2. ° 
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f b } are known to be subjected to the positivity constraint, i. e. > 

1 i 

b + B a 0 for all i. Several methods such as linear iteration, the 
i o 

method of steepest descents, and the method of conjugate gradients 
have been coded in Algol computer programs to solve the problem of 
nonlinear programming with linear constraints. 

3, 3. 1. Linear Iteration 

We can solve Eq. ( 2 . 52 ) directly by the method of linear itera- 
tion. From Eq. (2.52), we have 

t> 1+1 = A~ * Y +A _1 B(b 1 ) , 

where b* is the estimate of vector b in the i-th iteration. Matrices 
A and Y are not functions of b; hence they will not vary in the process 
of iteration. Only the B matrix is variable, and it is relatively simple 
to compute. 

The positivity constraint applies at each step in the iteration. 

A A 

Whenever b +B < 0 in the iteration, we set b. =. -B + e » where 
1 o 1 o 

e is a small positive number. The iteration process stops whenever 

the increment of b. is smaller than a predesignated value. The con- 

3 

vergence of the linear iteration scheme is rather slow, and it took 
approximately one hundred iterations to converge in typical cases. 
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3.3.2. The Method of Steepest Descents 

The linear iteration method converges rather slowly. In order 

to speed up the rate of convergence, we have tried second-order 

(6 7 8) 

methods such as the method of steepest descents. ' ’ This method 

has been described in section 2.5 to minimize F(b) of Eq. (2.47). 

In this second-order method we need to calculate the gradient vector 
V F and the Hessian matrix H of function F(b). The elements of the 
gradient vector and the Hessian matrix are 


e i = (v F) i 


_ _0F 
8b. 


N _ 


*-2 


i = l 


n. K. . 1 

nv. - “ ST K y 

o 1 


M 


2 

m= 1 J 


b , 
m m 


and 


c jk = H jk 


8 Z F 


8b. 9b. 

J k 


N 


v- n.K..K, 

* 2 -i-a lk 


} ' u jk 


i = 1 (I +V 
o i 


(3. 18) 


(3. 19) 


The values for e. and c.. are then substituted into Eq. (2.61) to 
i ij 

evaluate \ , and the new value for b is then obtained from Eq. (2. 59), 
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b. “ b , — X ^ * $ 3 b • « # > M • 

j j j 


i + 1 

The positivity constraint applies at every step; whenever , 

i + 1 

we set b. = -B + e . 

J o 

In this method, we have to update { e. } and [ c.. } in each 

1 i ij 

iteration step. This takes more effort than the linear iteration method, 
where we have to update only matrix B(b). The iteration stops when- 
ever 


M 


J | b. 1 +1 - b. 1 
j = 1 J J 


M 


£ 10 


-3 


I I b j l + 1 ! 

; - 1 J 


j = i 


The method of steepest descents converges quickly and start 

( 8 ) 

to oscillate when approaching the minimum point. 


3. 3. 3. The Method of Conjugate Gradients 

To remedy the disadvantage of the method of steepest descents, 
an attempt to use a more efficient computational method has been made. 
We have used the method of conjugate gradients. The method of con- 
jugate gradients has been well documented in the liter ature^^’ 
for simplicity we will only present the algorithm we have used. 

In the method of conjugate gradients, the major steps in the 


algorithm are: 
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1. At the initial point b* compute the gradient vector 
g, - VF(b^), and set the vector d, - -g, . 

^ I '''“"I -V 1 

2. For i = 2, 3, . . . , M +2 

(a) set b 1 =b 1 " 1 + X. . d 1 1 , where 

(b) X . , minimizes q. .(X) = F(b X ) ; 

l-l l-l ~ 

(c) compute g 1 = V F(b*) ; 

(d) unless i = M + 2 set 


,i i . ,i- 1 

d - -g + V d 

/V ^ 


where 


. i.T . i. i-l.T . i-1. 

Y = (g ) (g )/(g ) (g ) . 


and repeat (a). 

3. Replace b* by b^*^ and restart the iteration. 

We have applied this algorithm to minimize the functional F(b) 

>v 

with u -i - 0* By letting (j. = 0, we use no prior knowledge about 

j k j k 

the object. The elements { e. 1 } , { c * } of the gradient vector 

^ j / 2 \ 

g* = V F(b*) and Hessian matrix H* = < — — — ' ~ ^ i 

r, 


8 b i 8 b kfb= 


are shown in 


Eqs* (3. 18) and (3. 19) except with (j, = 0 C 

jk 

The above algorithm is for unconstrained minimization. Be- 
cause our functional is constrained, the algorithm has to be modified. 
The value of X in step 2(b) has been evaluated like the one in the 
method of steepest descents in Eq. (2. 61). We have applied the posi- 
tivity constraint in step 2(a), setting + e whenever 
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• i i-1 

b. 1 < -B and at the same time changing d. to (-B - b )/ X . The 

jo J ° J 

stopping criterion we have used is similar to the one used in the method 
of steepest descents. 

The positivity constraint presents difficulty for the method of 
conjugate gradients in our computer programming. The elements of 
conjugate gradient vector d have been changed whenever the posi- 
tivity constraint is violated; hence d 1 is no longer a conjugate gradi- 
ent vector. All the analysis for the convergence for this method is 
no longer necessarily valid. Up to now this method has not worked 

well. A better way to handle the positivity constraint such as penalty 

(13) 

function method probably will be needed. 

3.4. Computer Simulated Experiment Results 

In computer simulation, the one -dimensional object can be 
generated either by the first-order autoregressive process as de- 
scribed before or by just specifying any particular form. The parame- 
ters such as the uniform background level B q , the correlation coeffi- 
cient r, the signal mean deviation , and the constant asT should 
be given both in the processes of generation and restoration of the 
object. The parameters r and used in the restoration process 

are not necessarily equal to the values used in the generation of the 
object; usually the prior knowledge about them is not complete. The 
incoherent point-spread function used was that for either motion-blur 
image degradation or diffraction-limited distortion as described by 
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Eqs. (3.7) and (3. 12) respectively. 

The likelihood achieved by the maximum- likelihood estimate is 
higher than for the linear least-squares estimate for sure from the 
definition of maximum-likelihood estimation. If the mean error- 
squares of the maximum- likelihood estimate are smaller than for the 
linear estimate, then we can conclude that the maximum-likelihood 
estimate is superior to the linear estimate. 

In the nonlinear estimation scheme, we use both the linear 
iteration and the method of steepest descents to solve the set of non- 
linear equations (2.49). Prior knowledge of the positivity constraint 
on the object radiance function is used. We use the linear estimate 




k t k + ~ 

~ cc sT 




(3.20) 


as the starting point in solving the set of nonlinear equations. When 

2 

we assume the signal is much larger than the noise, i. e. , O *♦ 00 » 
the linear estimate in Eq. (3.21) will reduce to the inverse filtering 
estimate 


~ T -1 T / 1 \ 

b T = (K K) K I “ n - I I ) . (3.21) 

I ~ ~ ~ \ctsT ~ o ~ / 

The performance of different estimation schemes will be meas- 


ured by a percentage squared error from the formula 
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e 


M 

j = l 


(b. 


M 




X 100 


( 3 . 22 ) 


A 

where b, is the estimate of . We will compare the performances 
of the maximum- likelihood estimate, the linear least-squares estimate, 
and the inverse filtering estimate. For fair comparison the positivity 
constraint on the object radiance was also applied to linear least- 
squares and inver se -filter ing estimates. 

In a typical experiment, we used a set of five numbers (347, 
5412, 3298, 9547, 2486) as the initial number for the random number 
generator computer subroutine to generate random variables 
{ R* R 1 , Y. , R 1 , R 1 } . The object radiance function was de- 

-L b 1 J 

graded by motion blur with a blur of four cells, i. e. , D = 4 in Eq. 

(3. 7). The parameters used in the generation of the first object and 
data were (asT) = 10, = 3, a = 25 , r = 0. 6. The generated object 

and the measured data reduced to the scale of the object by dividing by 
(asT) are shown as the heavy line and the thin line in Fig. 3. 2 
respectively. Hard rectification has been applied at six object points. 
The measured data do not resemble the object at all. 

Because of the hard rectification, the object radiance function 
can no longer be described by a first-order autoregressive process 
having the covariance matrix as shown in Eq. (3.2). We used the same 
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values for the parameters r and 0^ in the restoration scheme as in 
the generation of the data. The result is shown in Fig. 3. 3. The 
signal-to-noise ratio is 5787 » 1, so the inverse filtering estimate 
is very similar to the linear least-squares estimate. 

In the restoration of the object radiance function, we have used 

values for parameters r* and that are not necessarily equal to the 

values r and a. we have used in the generation of the object radi- 
b 

ance function. From Fig. 3,4, we can see that different values used 
for parameters r' and affect the linear least-squares estimate 

more than the maximum- likelihood estimate. The maximum-likelihood 
estimator is a function of the measured data; as can be seen from 
Eq. (2.52), this estimator tends to be adaptive and hence depends less 
on the prior knowledge of the object radiance function than the linear 
least-squares estimator. This can be demonstrated by using a = <» 
in the restoration process; the restoration uses no prior information 
whatever about the object, but only the Poisson character of the dis- 
tribution of the data and the known point-spread function. The result 
is shown in Fig. 3.5. The percentage error-square for the maximum- 
likelihood estimate is 2. 69%, which is not too far from those results 
for which the prior knowledge was used. Another interesting point 
shown in Fig. 3.4 is that if the assumed correlation r is too large, 
say r = 0. 99, the result of maximum- likelihood estimation tends to 
go bad. This happens probably because if we assume that r is close 
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. Percentage squared error as function of assumed corre 
lation coefficient r ; curves are indexed by assumed 
standard deviation . Solid lines - linear estimate; 
dashed lines - maximum- likelihood estimate. 



Figure 3.5* Original object and restored values, x - 
likelihood estimate, o = linear estimate, 
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to 1 and actually it is not, the maximum- likelihood estimator will tend 
to closely correlate all the object points. In other words, the esti- 
mator has been very confined. In general, over a reasonable range, 
the bigger and the smaller r, the better the maximum- likelihood 

estimator performs. Under this situation, the maximum- likelihood 
estimator is not too confined and has the ability to adapt itself according 
to the measured data. 

The effect of asT is also examined. For an existing optical 
system, as is fixed, but the observation time T can be varied. The 
result is shown in Fig. 3. 6. Both the maximum- likelihood estimate 
and the linear least-rsquares estimate perform better for larger aST, 
since it is related to the s ignal-to-noise ratio, but the maximum- 
likelihood estimate is not as sensitive to the change of aST as the 
linear least-square estimate. 

Anotherinput object radiance function tested is shown in Fig. 
3.7(a). The parameters used in the generation of the object function 
are r = 0.2, CT = 15, and I = 3. The resulting restoration is shown 
in Fig. 3.7(b). Figures 3.8(a) and (b) show the percentage squared 
error as functions of various values of assumed correlation coefficient 
and standard deviation cj . The maximum- likelihood estimate is seen 
in this case to yield the closer fit to the true object radiance, the 
greater the assumed value of a and the smaller the assumed value 


of r . 
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Figure 3.7 (a). Original object (thick line) and its noisy version (thin 
line) as recorded after blurring by relative motion; 
(5431, 3214, 5438, 3654, 2495), 1 Q = 3, asT = 10, 
r = 0.2, a - 15. 
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Figure 3.7 (b). Original object and restored values. x = maximum- 
likelihood estimate, o = linear estimate, Q = . 
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Figure 3. 8 (a). 


Percentage squared error as functions of assumed 
object standard deviation a b ; curves are indexed 
by the assumed correlation coefficient r. Solid 
lines - linear estimate, dashed lines - maximum- 
likelihood estimate. 
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Figure 3. 8 (b). Percentage squared error as functions of assumed 
correlation coefficient r; curves are indexed by 
assumed standard deviation , Solid lines - 
linear estimate, dashed lines - maximum- likelihood 
estimate. 
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There is, however, considerable variability in the performance 
of both the linear and the maximum- likelihood estimators. We should 
like the estimator to perform well in most of the cases. We took the 
same blurred object, with the same set of mean values X ^ » generated 
a number of independent sets of data n. , and applied our procedure to 
each, with the resulting percentage errors shown in Table 1. The star 
indicates the trial shown in the previous figure. In ten of the twelve 
trials taking O = » gives a lower percentage error than when the 
restoration procedure uses the prior information about how the object 
is generated. In two trials out of twelve the linear least-squares filter- 
ing was more effective than the maximum-likelihood estimation. Most 
of the time, the maximum- likelihood estimate is considerably closer 
than the linear one. 

Table 2 shows a similar set of trials for a different object; 
there is much greater variability, but in all but one case the maximum- 
likelihood estimate was much closer than the linear one. Here in five 
out of twelve cases using the prior information gave a closer fit to the 
true object radiance. The starred case is shown in Figs. 3. 9(a) and 
3.9(b). 

Deterministic signals were also used. A typical object radi- 
ance function we used in the experiment has the form 

b^ ~ a{exp[-(i - 5)^/4] + exp[-(i - 9)^/4] + exp[-(i - 16)^/4 ] , 


(3.23) 
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Table 1 

Performances of both linear and maximum- likelihood 
estimation of an object based on twelve different sets 
of data. Parameters used to generate the object are 
(R^R^ = (5431, 3214), I q = 3, aST = 10, r = 0.2, O = 15. 


Squared Errors (percent) 


a = a , r = R o ® 

s' 



Linear 

Nonlinear 

Linear 

Nonlinea 

1 

5. 53 

1. 80 

8.81 

0.862 

2 

2. 13 

2. 46 

1. 02 

1.08 

3 

4. 90 

2. 93 

4.76 

2.21 

4 

5. 95 

4.28 

8.20 

2.44 

5 

7. 61 

5.20 

5. 17 

3.43 

6 

5. 83 

5.79 

3.43 

3.58 

7 

8.49 

4.23 

8.25 

3.68 

8 

6. 02 

5. 64 

3. 97 

3.69 

9 

5.49 

3. 67 

9.77 

4.76 

10 

9. 84 

6.22 

8.69 

5. 17 

11 

7. 84 

6. 36 

7.24 

5. 63 

12 

12. 59 

9. 39 

12. 36 

9.56 

Mean 

6. 85 

4. 83 

6.81 

3.84 

Standard 

Deviation 

2. 67 

2. 07 

3. 18 

2. 33 


Illustrated in Fig. 


3.7. 
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Table 2 

Performances of both linear and maximum- likelihood esti- 
mation of an object based on twelve different sets of data. 
Parameters used to generate the object are (RpR^) = 
(842, 3248), I* = 5, aST = 10, r = 0. 3, a = 15. 




Squared Errors 

(percent) 



a ~ a , 

s 

r = R 

(J -> 00 



Linear 

Nonlinear 

Linear 

Nonlinear 

1 

7.91 

6.70 

7. 74 

5.42 

2 

9.69 

7. 07 

8. 11 

5.70 

3 

14.93 

8. 84 

11.82 

6.59 

4 

13.57 

7.23 

15.79 

6.90 

5 

9. 16 

5. 10 

11. 55 

7.51 

6 

13. 94 

8. 12 

18.44 

9.73 

7 

19.62 

11.75 

19. 91 

* 

10. 17 

8 

9.51 

6.86 

9. 94 

10.74 

9 

29. 62 

8.62 

47. 06 

10.99 

10 

16. 04 

14.55 

16.81 

14.27 

11 

23. 18 

18. 07 

24.70 

17. 02 

12 

22. 19 

12. 78 

25. 73 

17. 05 

Mean 

15.78 

9. 64 

18. 13 

10. 17 

Standard 

Deviation 

6. 70 

3. 85 

10. 90 

4. 11 


$ 

Illustrated in Fig. 3.9. 




Figure 3. 9 (a). Original object (thick line) and its noisy version (thin line) as re 
corded after blurring by relative motion; (842, 3248, 854, 

279, 657). 1 Q = 5, asT =10, r = 0. 3, a = 15. 


sO 

CO 
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where a is a positive number which corresponds roughly to the stand- 
ard deviation of the object. The object radiance function consists of 
three peaks of height a; two peaks are close together and one peak is 
more isolated. The width of each peak is about four. For a = 90, 
the object and the smeared image are shown in Fig. 3. 10. The object 
has been smeared by a motion blur of four cells. It is clear from 
Fig. 3. 10 that two peaks have been smeared into a broad peak. The 
distance between those two peaks is four cells, which is the amount of 
the motion blur. 

A set of data using parameters asT = 3 and 1=5 was gener- 
ated for the Gaussian- shaped object with a = 90. A typical restoration 

result using parameters r' = 0.5 and o 7 = 45 is shown in Fig. 3. II. 

s 

Two peaks are well resolved for both linear and maximum- likelihood 
estimate. The performances of both estimators as functions of (aST) 
for the same object are shown in Fig. 3, 12. The effect of using differ- 
ent assumed values of r 7 and is shown in Figs. 3. 13 and 3. 14. 

Even though for some values of r 7 and , the linear estimate per- 
forms better than the maximum- likelihood estimate, in most of the 
cases the maximum- likelihood estimate performs better than linear 
estimate. For deterministic signals, all the statistical assumptions 
about the object are not valid any more; however, we can still have 
good results. This demonstrates the usefulness of our estimator. 

Image degradation by a diffraction-limited optical system has 







PERCENTAGE 
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Figure 3. 14. Percentage squared error as function of assumed 
correlation coefficient r . CtsT ~ 3, - 5, 

a = 45, a = 90. Solid line - linear estimate; 
dashed line - maximum likelihood estimate. 
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also been considered. Figure 3. 15 shows the result of the restoration 
of a stochastic signal. The original object has been discussed in sec- 
tion 3. 1 and shown in Fig. 3.1. The passband width W q in Eq. (3. 13) 
has been set to 1.2. In this case, the object spectrum outside the 
passband of the optical system cannot be recovered; hence the small 
sharp peak between two big peaks is not restored. Figure 3. 16 shows 
the performances of both linear and maximum- likelihood estimators 
as functions of various assumed values for r and . For = 120, 

r = 0. 8, the percentage squared errors are 11.72 and 32.60 for the 
nonlinear and linear estimators respectively; and for smaller r the 
performance of the maximum- likelihood estimator deteriorates some- 
how, but the percentage squared error exceeds 100%. Owing to the 
adaptive characteristic of the maximum- likelihood estimate, its per- 
formance is rather insensitive to the assumed value of the parameters 

used. We have used a ' = 00 in our restoration scheme. The percent- 

b 

age squared errors for the nonlinear estimator are 24. 52% and 4880% 
for the linear estimator. 24. 52% is about the best the linear estimator 
can do as we can see from Fig. 3. 16. The nonlinear estimator always 
performs better than the linear estimator in this case. 

Figures 3. 17 and 3. 18 show the performance of the restoration 
of a deterministic Gaussian- shaped object. The object has been de- 
graded by a diffraction-limited optical system with W q = 0. 8. It is 
clear that maximum- likelihood estimator performs better than the 


linear estimator. 



2 4 6 8 10 12 14 16 18 ~ 20 


Figure 3. 15. Original object and restored values, a = maximum likelihood 
estimate, o = linear estimate, r = 0.2, CT = 20, W Q - 1.2. 
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Figure 3. 16. Percentage squared error as functions of assumed 
correlation coefficient r . Curves are indexed 
by the assumed object standard deviation <7 . 

Solid lines - linear estimate; dashed lines - 
maximum- likelihood estimate. 
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Percentage squared error of the restoration of a 
Gaussian- shaped object as function of CtsT . 

Solid line - linear estimate; dashed line - maximum 
likelihood estimate; O- 15, I Q = 10, r =0.5, 
w = n a 


PERCENTAGE ERROR 



In summary, we have developed a method of restoring images 
that takes account of the Poisson character of the recorded image data, 
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instead of making the assumption of Gaussian-distributed data that 
underlies the usual linear least-squares procedure. The significance 
of the new restoration procedure has been demonstrated by computer 
simulations. In the great majority of trials the new method produced 
a much better fit to the true object radiance than the conventional 


linear procedure. 
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Part II. Optical Signal Processing by Shearing- 
Interferometric Methods 

In part II, we are concerned with signal processing using 
coherent optical systems. 

Communication theory has been successfully applied to optical 
signal processing, particularly in the areas of holography and spatial 
filtering. Many optical phenomena can be easily explained using com- 
munication theory. Carr ier -frequency photography is an amplitude- 
modulation (AM) process; computer generated holography is just a 
pulse-width (PWM) and pulse-position modulation (PPM) scheme. 
Diffraction by a grating is a harmonics generation process, and the 
moire effect can be understood by modulation theory. Modulation, 
demodulation, heterodyning, generation of harmonics, and filtering 
have been integrated into communication theory. In optics, however, 
there are as yet few connections between grating diffraction, the moire 
effect, carrier-frequency photography, holography, and interferometry. 
Here we want to integrate these processes and apply them to optical 
signal processing. 

Following an introduction and a brief review of the pertinent 
details of grating diffraction, Ronchi's interferometer, the moire 
effect, carrier-frequency photography, and achromatic holography, 
two types of interferometers are presented. In these interferometers, 
diffraction gratings are used at either the spatial-frequency plane or 
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the object plane of a telecentric image-forming system. In the grating 
interferometer, gratings are placed at the spatial-frequency plane as 
beam splitters, and in the carrier -frequency interferometer gratings 
are placed at the object planes of two cascade telecentric image- 
forming systems as modulation and demodulation masks. These inter -_ 
ferometers can produce variable shearing in both lateral and radial 
directions. The operation of the grating interferometer is 
in real time; the operation of the carrier -frequency interferometer can 
be either in real time or in two steps. 

Optical testing, phase-object visualization, spatial differentia- 
tion, complex amplitude addition and subtraction, image contrast re- 
versal, image complementary color reversal, two-dimensional 
photometry, fabrication of filters, pure spatial sinusoidal wave genera- 
tion, and image multiplexing can all be achieved by using these new 
interferometers. The methods employing the new interferometers have 
some definite advantages in terms of simplicity, flexibility, reliability, 
light efficiency, or real-time operation in comparison with existing 
methods . 

Optics and electronics are means of information processing 
which have many common features. Communication theory is well 
developed in electronic signal processing. Recently it has been suc- 
cessfully adapted by optical physicists, especially in the development 
of holography and spatial filtering. Many optical phenomena can be 
easily understood from the point of view of communication theory. 
Modulation, demodulation, heterodyning, harmonics generation, and 



110 


filtering have been used successfully in electronic signal processing; 
we can find their counterparts in optics. In optics, modulation can be 
achieved by carrier -frequency photography; generation of harmonics 
by grating diffraction; demodulation and heterodyning by the moire 
effect; and filtering by spatial filters. Carrier-frequency photography, 
diffraction by gratings, and the moire effect all have their unique func- 
tions in optical signal processing. Two new types of interferometers 
have been developed from the combined use of those concepts. 

Diffraction by gratings and the moire effect have been used ex- 
tensively in optics. Gratings have been used as beam splitters. 

Ronchi invented a shearing interferometer using a grating as a beam 
divider, which has been used extensively for testing optical components 
such as lenses. This technique is relatively simple and inexpensive. 
But the use of just one grating results in many drawbacks, such as 
(a) the fixed amount of shear, (b) the dependence of the fringe contrast 
on the degree of equality in light flux between different diffraction 
orders, and (c) the overlap of many diffraction orders. 

The drawbacks of Ronchi 's interferometer can be overcome by 
the use of two superimposed gratings instead of just one grating as a 
beam divider. The amount of shear in Ronchi's interferometer is pro- 
portional to the spatial frequency of the diffraction grating used. If we 
can vary the spatial frequency, we can have variable shear. This can 
be done by the moire effect of the two superimposed gratings. When 



Ill 


two gratings are superimposed, owing to the moire effect a variable 
spatial frequency is produced by rotating the angle between their axes. 
An interferometer, called the "grating interferometer, " has been 
developed from two superimposed Ronchi rulings used as beam split- 
ters. All the drawbacks in Ronchi's interferometer are overcome. 

Carrier-frequency photography has some similarity to holog- 
raphy. In particular, both are two-step processes, and diffraction 
methods are generally used in the retrieval stages. Holography has 
been used very successfully in interferometry. Holographic inter- 
ferometry can be operated either in real time or in two stages. 

A hologram's ability to record both the amplitude and phase of 
a wavefront is the reason for its success in interferometry. Carrier- 
frequency photography can register the amplitude as well as the phase 
variation of objects. When a grating is imaged through a complex ob- 
ject, the resultant image is an amplitude -modulated version of the 
grating with the slits shifted according to the phase variation of the 
object. This kind of encoding scheme, involving both amplitude - 

modulation and pulse-position modulation, called the "detour phase" by 

(2) 

A. Lohmann, is the essence of computer -generated holography. It 
is very difficult to retrieve the information about the phase variation of 
the object directly from the distorted grating image. But if we "beat" 
the distorted grating image with a regular grating (or master grating), 
moire fringes appear which will reveal the phase variation. A new 
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type of interferometer, based on carrier-frequency photography and 
the moire effect, has been conceived. It is called the carrier -frequency 
interferometer* It can be operated in real time or in two stages. The 
system is achromatic and white light can be used. 

In chapter 4, we briefly review diffraction by gratings, the 
moire effect, carrier -frequency photography, and achromatic holography. 
This review will serve as the foundation of the following chapters. 

The grating interferometer is the subject of chapter 5, and the 
carrier-frequency interferometer is described in chapter 6. Many 
kinds of masks can be used as modulation carriers in the carrier- 
frequency interferometer. We use linear gratings, checkerboard 
gratings, and circular gratings. 

The applications of these two types of interferometers in optical 
information processing are presented in both chapters 5 and 6 for each 
interferometer. New methods for optical testing, phase object visual- 
ization, image polarity reversal and color complementarity, two- 
dimensional photometry, spatial differentiation, complex amplitude 
addition and subtraction, image multiplexing, pure spatial sinusoidal 
wave generation, fabrication of filters, focusing, and alignment are 
presented. Some of the applications can be used in optical computing, 
general image processing, and two-dimensional photometry. 



Chapter 4 


Review of Diffraction by Gratings, the Ronchi Interferometer, 
the Moire Effect, Carrier -Frequency Photography, 
and Achromatic Holography 

Gratings have long been used as beam splitters, particularly in 
the Ronchi interferometer. The moire effect has been applied to 
measurement processes, carrier-frequency photography has been ap- 
plied to image multiplexing, and holography has been applied in inter- 
ferometry. They all have their own unique applications to optical 
information processing. However, few efforts have been made to inte- 
grate all these techniques. It is the goal of this part of the thesis to 
integrate them for use in optical information processing. In this 
chapter we will review them briefly. 

4-1. Gratings, Diffraction by Gratings , and Ronchi Interferometer 
One type of grating can be represented as a periodic binary 
mask. A grating with slits parallel to the y-axis can be represented 
as 


j 1, 1 (x - a)/d - N| £ q /2 

M(x) = < , N = 0, 1, 2, . . . (4. 1) 

( 0, otherwise 

where d is the period of the grating, a is the shift from the origin of 
symmetry, and q is the ratio of the size of the openings (of unit trans- 
mittance) to the period of the grating. For the Ronchi ruling, q - 1/2, 
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so that the size of the opening is equal to half the period of the grating. 

M(x) is a periodic function, which can be expressed by means 
of Fourier series, ' 

M(x) ^A n e Z " i " t/d , (4.2) 

n 

where = q sinc(nq) exp(2 rr i na/d) , 
sinc(x) - sin(TT x)/{rr x) 



(here * denotes the complex conjugate). The binary grating M(x) 
and its Fourier series spectrum for the case when q = 1/2, are shown 
in Fig. 4. 1. No even harmonics exist. 

If we rotate the grating by an angle 9 about the origin, its 
transmittance function becomes 

M(*. y) = ^ A exp [2rri n(x cos 0 + y sin0 )/d] . (4.3) 

n n 

The grating has the frequency cos 9 /d in the x-coordinate and sin9 /d 
in the y-coordinate. 

The periodic binary mask M(x) can be generalized to include 
quasi-periodic masks, e. g. the circular grating. 


From now on, ^ means summation over n from - 00 to +® . 
n 



M(x) 




Figure 4. 1. A Ronchi ruling and its Fourier spectrum . 

(a) The transmittance function of a Ronchi ruling. 

(b) The Fourier coefficients of the Ronchi ruling; 
no even harmonics exist. 
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( 1 , | r/a - N| * q/2 

C(r) = j N = 0, 1,2, ... (4.4) 

l 0 , otherwise 

When a monochromatic plane wave in z-direction is used as 
incident wave, the wavefronts diffracted by a transparent one- 
dimensional plane grating have their wave vectors in the plane of both 
the original wave vector and the direction perpendicular to the direc- 
tion of the slits (or in the direction of grating orientation). The grating 
splits a wavefront incident on it into many diffraction orders, accord- 
ing to the relation ml/d = sin 9 m » where m is the diffraction order, 

X the wavelength, 8^ the angle between the zeroth order and the 
m-th order. The amplitude of m-th diffraction order is . 

When a transparent one -dimensional plane grating is placed in 
the spatial-frequency plane of a telecentric image -forming system as 
shown in Fig. 4.2, the grating acts as an amplitude -modulation mask 
in the spatial-frequency plane. The grating will generate harmonics of 
the object. By Fourier analysis, the field propagated through the 
grating will split into many orders in the image plane. Ronchi invented 
an interferometer using a grating as beam splitter. Figure 4. 2 is a 
typical setup. A monochromatic plane wave illuminates the object. 

The field just before the grating, owing to the Fourier transforming 
property of lens Lj , is 

U (x, y; 2 fj - 0) = U Q (x/ X f j , yf X f x ) , 


(4.5) 



SOURCE 


Figure 4.2. A coherent optical signal processor (or telecentric image-forming system). 

OBJ is the plane where the input object will be located, F is the filter plane 
where the spatial filter is located, and IMG is the image plane, and {2 
are the focal length of lenses Lj and respectively. 
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where 
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U Q (x, y) = JJ U Q (x / 1 y') exp[-2rri(xx' + yy ' )] dx' dy' . (4. 6) 

We have omitted an unimportant multiplicative constant. Applying the 
Kirchhoff boundary condition to the field as it passes through the grating, 
we have 

U(x, y; 2 f x + 0) = U(x, y; 2 ^ - 0) G(x, y) 

= U Q {x/Xf 1 , y/X fj) ^ A n exp(2rr i nx/d), (4.7) 

n 

where G{x, y) = ^ A exp(2rri nx/d) is the grating transmittance func- 
n 

tion. The field propagates through lens and reaches the image 

plane. Owing to the Fourier transforming property of lens and the 

convolution theorem of Fourier transformation, the field at the image 
plane is 

U(x, y; 2 ^ + 2 f^ = U(x, y; 2 ^ + ) 

= ^A n U 0 (-f 1 f' 2 1 3£-nXf 1 d‘ 1 l y). (4.8) 

n 

The object field has been reversed, magnified f^/f^ ^ mes > an d dif- 
fracted into many orders with amplitude A^ for the m-th diffraction 
order, and shifted over a distance . 

If the amount of shear, Xfj/d, is smaller than the size of the 
object, fields coming from different diffraction orders will overlap. As 
all the fields from different diffraction orders are coherent, interference 
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fringes will appear in the overlapped area. This interferometer has 
been used extensively in optical testing. 

4.2. The Moire Effect 

The moire effect was described by Lord Rayleigh^ ^ as long 
ago as 1874. When two high-spatial-frequency gratings are superim- 
posed, by proper arrangement, a pattern with low-spatial-frequency 
variation appears. 

Two linear gratings T^x, y) and T^x, y) have transmittance 
functions 

T x (x, y) = ^ A n exp [2TTin( Vl x + Mjy) ] , 
n 

and (4. 9) 

T 2 (x, y) = ^ B m exp [2ni m(v 2 x + (i 2 y ) ] , 

m 

where (Vj.Hj) and (v 2 » are fundamental spatial-frequencies in 
{x, y) of gratings T^ and T 2 respectively. When they are superim- 
posed, the resulting transmittance is 

T (x, y ) = TjCx, y) T 2 <x, y) 

L B m eXp ^ rri t^ nV i + mv 2^ x + ^ i + 10 ) 

m n 

Examining the Fourier spectrum of the transmittance functions 
of Tj , T 2 , and T, we see that T^ has the fundamental frequency 
(v,|j) = (Vj.Mj). T 2 has (v 2 , M- 2 ^’ ant ^ resu l t l n g transmittance T 
has (Vj-V2 * |i i“M 2 ^» which is the beat frequency of Tj and T 2 . 
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4.2. 1. Moire of Two Linear Gratings of Different Period, 
Parallel Orientation 

Two gratings T. and T_ both have their slits parallel to 

X c* 

y-axis, and have spatial frequencies Vj and respectively. The 

total transmittance of these two gratings superimposed, by eq. {4. 10)), 
is 

2TTi(nv, + mv 9 )x 

e 1 L (4.11) 

m 

m n 

When both T. and T are Ronchi rulings, q = 1/2, no even harmonics 

X L* 

exist. For the case that w * an< ^ PI = 1/2, there are only two 
diffraction orders (n=l, m= 0) and (n = 0, m =1) close to v = jT (v^ + V^), 
as shown in Fig. 4. 3. 


T(x, y) B 


4.2.2. Moire of Two Linear Gratings with Same Frequency, 
Different Orientations 

When the two gratings T^ and T^ have the same period d, 
but their orientations are different, one grating having its slits rotated 
through an angle 9 from the x-axis, and the other having its slits 
rotated through an angle -0 from the x-axis, the total transmittance is 


-H? 


T(x # y) = < / A exp[2tT i n(x sin 0 + y cos 0 )/d] 

1 / - n 


(2 


m 


exp[2iT i m(-x sin 0 + y cos 0 )/d] 


I 

n 


2 


m 


A B exp tetri (xv 

n m * L nm 


+ yu ) 1 

7l nm J 


( 4 . 12 ) 
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Figure 4. 3. The Fourier spectrum of two superimposed Ronchi 
gratings of different periods, parallel orientation. 
There are only two diffraction orders (1, 0) and 
(0, 1) that are nonzero and satisfy m + n - 1. 


Figure 4.4. The Miller's indices. For Ronchi rulings there are 
only two diffraction orders that are nonzero and 
satisfy m + n = 1 and -1. 
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where 

V = (n - m) d 1 sin 8 

nm 

u = (n + m) d 1 cos 0 . 
nm 

The Fourier spectrum of T(x, y) can be visualized in the (v,|j) 
domain with the aid of "Miller's indices" (n, m). Figure 4.4 shows the 
location of different indices in the (v , ) plane, for the case that both 

T and T are Ronchi rulings. The amplitude strength of each dif- 

i. 4 

fraction order (n, m) is A B . All the diffraction orders with 

n m 

n 4- m = p are located at jj = pd cos 0 , with a distance 2d sin 0 
between the neighboring diffraction orders. If we vary Q , the distance 
2d” 1 sin0 is varied. The distance 2d" 1 sin 0 is the beat frequency 
produced by T^ and T superimposed. When both T^ and T^ are 
Ronchi rulings, no even harmonics exist, i. e. , A = 0, for n or 
m an even number. There are only two diffraction orders, (0. 1) and 
(1. 0), with m + n = 1 , as shown in Fig. 4. 4. 

4.2. 3. Moire of One Regular Grating and One Distorted 
Grating 

When a regular grating G^(x, y) = a +b cos (vx) and a distorted 
version of that grating G^fx^y) = a + b cos [vx + 9 (x, y)] are super- 
imposed, the total transmittance is 

T p (x, y) = Gj (x, y ) G 2 (x, y ) 

2 

= ^ b^cos 0 +2 ab cos|| jcos (vx+^J + cos (2 v x + 0) , 

(4. 13) 
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where 0 (x, y) is related to the distortion of position of the grating slits. 
The information which we wish to observe, 0 (x, y), can be extracted by 
spatial filtering if we let pass only the first diffraction order in 
Eq. (4. 13). 

If we add the transmittance functions of and together, 

the total transmittance is 

T (x, y) = 2a + ^ cos ^ v x + |- j cos . (4.14) 

The addition of transmittances can be achieved optically by a double - 
exposure technique. The information about the distortion function 
appears as the envelope of the carrier wave cosfvx + 0/2). 

4.3. Carrier -Frequency Photography 

In communication, modulation theory has been used in connec- 
tion with the transmission of information. Several signals can be 
transmitted over a single cable at the same time, by either time multi- 
plexing or frequency multiplexing. In optics, we can do essentially the 
same thing for multiplex information storage, by which we mean that 
more than one image is stored in a single film frame. By analogy, 
optical image multiplexing can be achieved by either spatial multiplexing 
or spatial-frequency multiplexing. Spatial multiplexing means each 
single image occupies a certain part of the film frame, which can be 
done by either a microfilm technique or image interlacing. In optical 
spatial-frequency multiplexing, each single image is carried by carriers 
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of different spatial frequencies. As a different orientation of the same 

grating yields different spatial frequencies, we have a theta-modulation 

, ( 2 ) 

system. 

Holography can record both phase and amplitude of the object 
wavefront and has been successfully applied in interferometry. Carrier- 
frequency photography has some similarities to holography, but the 
recording procedure is different. We will apply carrier -frequency 
photography to interferometry. 

The basic recording and reconstruction setups for carrier - 
frequency photography are shown in Fig. 4. 5. In recording, the object 
is imaged onto a film, and the modulation mask (a Ronchi ruling) can 
be placed in either the film plane or the object plane. We can build a 
camera with a grating in contact with the film plane. In reconstruction, 
coherent light is used to illuminate the film, and a binary filter is 
placed in the filter plane. This binary filter passes just one of the 
diffraction orders. 

If we want to store more than one signal, e.g. four signals of 
intensity distribution functions I {x, y), I {x, y), I (x, y) and I (x, y), in 

J- j *x 

(3) 

a single film frame, we can use theta modulation. Four exposures 
are made, each exposure for a different object. Between each exposure, 
the grating is rotated through an angle 0 = 45 ° , hence the designation 
"0 - modulation. " The total exposure is 
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OBJ G + FILM 



(c) 

Figure 4.5. The recording and reconstruction setup for carrier- 
frequency photography, (a) and (b) are for recording 
and (c) is for reconstruction. The filter in (c) passes 
only one diffraction order. 
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I(x, y) = I^(x, y) G(x, y; 0) + 1 ^ (x, y) G(x, y; tt/4) + I^( x » y) G(x, y; tt/2 ) 

+ I^(x, y) G(x, y; 3rr/4) , (4.15) 

where G(x, y; 0 ) is the grating transmittance function rotated through 
an angle 9 . When the film is developed in the linear portion of the 
characteristic T-E curve, the amplitude transmittance of the film will 
be proportional to I(x, y). 

In reconstruction, the developed film is placed in the object 
plane of a telecentric image forming system as shown in Fig. 4.5(c). 

The spectrum of I(x, y) in the filter plane will be similar to the one 
shown in Fig. 4. 6. The spectra of different objects occupy lines of 
different angles 9 , each split into many diffraction orders. Each sig- 
nal can be retrieved separately, or the signals can be combined by using 
a binary filter to pass only the spectra of the images to be displayed. 

This technique has been successfully applied to color photography using 

(3) 

only black and white film. 

Image multiplexing can be achieved by spatial interlacing as well 
as by theta modulation. In spatial interlacing, the grating is displaced 
a distance qd between each exposure; q as defined before is the ratio 
of the slit width to the period of the grating. We can record q 1 signals 
and retrieve them separately. For a Ronchi ruling, q = 1/2, we can 
record two signals in one film frame. The total transmittance is 



Figure 4. 6. The spectrum of a theta-modulated carrier- 
frequency photograph. 


Figure 4. 7. A double-exposure carrier-frequency modulated 
photograph. 
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I (x, y) = Ij(x, y) G(x, y) + G ^ x " d ^ 2 ’ ^ 

= ^ A [I.(x, y) + (- 1 ) n I_(x,y)] exp{2TT i nx/d) . (4.16) 

n l £ 

n 

I(x, y) is shown in Fig. 4. 7. In reconstruction, we can put the modu- 
lation mask in contact with the developed film to pass only either 

I (x.y) G(x, y) or I (x, y) G(x - d/2, y), and retrieve them individually. 
L £ 

Or if we use no mask, and if the filter F passes only the first diffrac- 
tion order, | I (x, y) - I 2 (x,y)| 2 results. We create bipolarity by 
spatial interlacing. This technique has been used in image subtrac- 
tion^ and image polarity reversal. ^ 

4. 4. Achromatic Holography 

Because of the annoying twin image, holography was not very 
popular after Gabor invented it in 1948. The renaissance of holography 
came after Leith and Upanieks showed that the twin image can be 
eliminated by the use of an offset reference wave. The typical setup 
is shown in Fig. 4. 8. A point source of illumination is collimated by 
the lens L. The beam splitter BS divides the beam into two parts, 
with one part of the beam striking the object, which is taken to be a 
transparency with a complex transmittance 0(x, y). The second por- 
tion of the plane wave strikes a mirror and reflects into the film with 
an angle 0 with respect to the optical axis of the lens L . 
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Figure 4.8. A recording setup for holography. BS is beam 
splitter, HO IX) is hologram plane, M is mirror. 


G +OBJ F H 



Figure 4. 9. A recording setup for achromatic holography. The 
filter F consists of a pinhole at the center of the 
zeroth diffraction order and a binary filter at the 
first diffraction order. 
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The amplitude distribution across the film can be written as 


U(*,y) = R e 2T ' ia ' r + a<x,y)e +i * (X - y) 


where R is the amplitude strength of the reference wave, a(x, y) and 
$(x, y) are the amplitude and phase of the wave originating from the 
object, and a = sin 0/X is the spatial frequency of the reference wave. 
The intensity distribution across the film is 


E(x, y) = | U(x, y)| 2 

2 2 

= R + | a(x, y) | + 2 R a(x, y) cos [2ttoc y - 0(x, y) ] . 

Since a (and possibly a{x, y) and 0{x, y)) is dependent on 
wavelength X , each wavelength X produces its own elementary fringe 
system, and the systems of different wavelengths are shifted with 
respect to each other. Monochromatic light is necessary. 

If the interference fringes on the hologram are achromatic, a 
broad-band source such as natural white light emitted from a xenon-arc 
lamp can be used as illuminating source. To produce an achromatic 
hologram, it is sufficient that (1) the interference fringes on the holo- 
gram be achromatic, and (2) the diffraction pattern from the object be 
fairly independent of wavelength X in the plane of the hologram. ^ 
These requirements can be achieved by the optical setup shown in 
Fig. 4. 9. 

The object Uq(x, y) Is modulated by a grating G. By the 
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Fourier transforming property of lens , the field at the filter 
plane F, for each wavelength is 

U(x,y; X; 2f-0) = A r U q (x - n X f /d, y ) , (4.18) 

n 

where f is the focal length of both lenses and J-^ • If the illumi- 

nating source is polychromatic, each diffraction order except the 
zeroth order will exhibit color dispersion. We can use a binary filter 
F to pass the +1 diffraction order as the object wave, and a small 
pinhole in the center of the zeroth diffraction order to create the refer- 
ence wave. The field after the filter plane is 

U(x, y; 2f +0) = R8(x) + A 1 U Q (x- \fd~\ y) , (4.19) 

where R is the strength of the reference wave. The object waves pro- 
duced by light of different wavelengths are shifted with respect to each 
other in the filter plane F, but they will be recombined if the wave 
propagates to the hologram plane H as shown in Fig. 4. 9. The total 
field in the hologram plane H is 

U(x, y; X ; 4f) = R + Aj U Q (-x, -y) exp(2TTix/d) , (4.20) 

which is X -independent. Thus an achromatic hologram is recorded. 

White light can be used in holography. We will use the achro- 
matic holographic technique in a new type of interferometry, which can 
be used in several ways, especially for real-time image contrast re- 
versal. Positive black- and-white images can be obtained directly from 
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color negatives in real time without complicated photo -chemical and 
electronic techniques. 


\ 
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Chapter 5 


Variable Shearing Interferometry Based on the Moire 
Effect of Gratings - The Grating Interferometer 

5.1. Introduction 

Shearing interferometry has been used to test instruments and 
optical components. In this type of interferometry, the wavefront under 
test is duplicated and displaced. The two resulting wavefronts are made 
to interfere with each other, yielding a comparison of the wavefront 
with a shifted image of itself. These interferometers typically used a 
single beam-dividing element that is simpler than other interferome- 
ters. ^ There are two ways of obtaining the shifted wavefronts. The 
first method is to have the shift sufficiently large so that the displace- 
ment between the duplicated wavefronts is larger than the dimension of 
the wavefront under test. This is called "total shear, " since the 
amount of shift is larger than the size of the object wavefront under 
test. In the second method the shift is kept smaller than the width of 
the typical details of the object under observation; this is known as 
"differential shear, " owing to the similarity to differentiation. It is 
desirable to have an interferometer which can vary the amount of the 
shift. 

The interferometer invented by Ronchi uses a grating as the 
beam divider. This interferometer is simpler and cheaper than most 
other methods and has been used extensively. However, as mentioned 
before, the use of only one grating causes some inconveniences. These 
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can be overcome by using two gratings as a beam splitter instead of 

just one. Lohmann and Bryngdahl have described an interferometer 

using two diffraction gratings in the spatial frequency plane of two 

( 2 ) 

cascaded telecentric image-forming systems. 


5.2. Principle of the Interferometer 

We are proposing a new shearing interferometer based on the 
moire effect of two superimposed gratings. The optical arrangement 
of the grating interferometer is schematically illustrated in Fig. 5. 1. 

A plane wave falls onto the object plane OB J, which is imaged 
via telecentric image-forming lenses and onto the image 

plane IMG. The two superimposed gratings Gj and G^ are placed 
in the spatial-frequency plane. Both gratings Gj and G^ are linear 
and have the same period d; one has its slits rotated through an angle 9 
with respect to the x-axis; the other has its slits rotated through an 
angle - 0 . The transmittance of the superimposed gratings, as 
described by Eq. (4. 12), is 

T <*.Y>= 22 A n B m eXp P T ' i<XV nm +yw nm )]l <5> l) 

m n 


where 

V = (n - m) d * sin 0 
nm 

u = (n +m) d * cos 0 , 
nm 

and A and B are the Fourier coefficients of gratings G. and G_ 
n m 1 2 




Figure 5. 1. The grating interferometer. 
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respectively. The wavefront diffracted by this composite grating will 
split into many diffraction orders with the centers of each diffraction 
order distributed like the "Miller's indices, " as shown in Fig. 4.4. 
Variation of the angle of rotation 8 changes the distance between these 
centers, and hence the amount of shift between the different diffraction _ 
orders. 

The object U Q (x, y) is illuminated by a collimated plane wave. 

The field in the plane just preceding the grating, by the Fourier trans- 

(3) 

forming property of lens , is 


U (x, y; 2 f j - 0) = U Q (x/ \ fj , y/X fj) (5. 2) 

where f ^ is the focal length of lens . By the application of the 
Kirchhoff boundary condition, the field just behind the grating is 

U(x,y; 2^+0) =£ £ A n B mV x/X f l’ y/X V exppni,X V + nm» * 

n m 

(5. 3) 

This field will be once again Fourier transformed when it propagates 
through lens and arrives at the image plane. By the use of the 

shift property of Fourier transformations, the field in the image plane 
can be written as 


U(x,y; Zf 1 + 2f 2 ) = £ £ 

n m 


1 

AB U I - j — x - v 
n m 0 l 1-. nm 


, , 1 

x f r ‘ t y _|a 
2 


X£, 

nm 1 


(5.4) 
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where f is the focal length of lens L . Equation (5.4) can be written 
as 


U(x,y; Zi 1 +2f 2 ) = ^ U(x,y; p) , 


(5.5) 


and 


U(x,y; p) =V A B U.f-M^x-nA , -M^y-pA), 

y ^ n+p p, n o x y 

n ~T 2 


(5.6) 


M = f z /£ 1 , 

A = Xf, d" 1 

X 1 

A = Xf, d" 1 

y i 


sin 0 


cos 0 


and only A. ’s and B.'s indexed by integers are meaningful; when p 

is an odd number, n must be an odd number, in order that A , B 

njp p-n 

2 2 

be meaningful. For diffraction orders located at the same y- 
coordinate, the distance between each neighboring order is 2M A = 

X 

2 X f 2 d * sin 9 . All the diffraction orders U(x, y; p) appear at the 
same time in the image plane with their y- coordinate shifted by -MpA . 

The object field in the image plane has been inverted and magni- 
fied by a factor M = f^/fj , and the center of the diffraction order 

(n, m) is at x = (n - m) M A , y = (n +m)M A . If the dimension 6 

x y y 

of the object in the y-direction is smaller than A , all the diffraction 
orders with m + n = p will appear separately from those with 
m + n = p' , p' ^ p. We have achieved total shift in the y-direction. 
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If we use Ronchi rulings, q - 1/2, then no even harmonics exist, i. e. , 

A = b =0, for n even. For diffraction order p = n + m = 1, {n, m) 
n n 

can have values (1,0), (0,1), (2,-1), (-1,2). Only diffraction orders 
(1, 0) and (0, 1) do not have even harmonics, hence can be nonzero, and 

U(x,y; p = 1) = A Q B 1 U 0 (-M" 1 x +A x , -M _1 y - A y ) 

+ A 1 B () U 0 (-M" 1 x - A x , -M _1 y - A y ) . (5.7) 

The U(x, y; p = 1) diffraction order is centered at y = -M . The 
distance between the (0, 1) and (1, 0) diffraction orders is 2 M A^ , 
which varies with the angle 0 . The new interferometer is thus a 
/ variable shearing interferometer with only two diffraction orders over- 
lapping. If gratings and G 2 are identical , | A n | = 1 | , 

hence 1 A^BjJ = | A^qI ’ t ^ ie two di ff raction orders (0, 1) and (1,0) 
have the same strength. The phase difference between A^ and 
can be adjusted by shifting either G 1 or G 2 in the direction perpendic- 
ular to its slits. All the inconveniences of the Ronchi interferometer 
are overcome. 

The total output image field is represented by Eq. (5.5). The 
p = 1 diffraction order represents a two-beam interferometer; so does 
the p = -1 diffraction order. The p = 1 diffraction order is centered 
at y = MA (compare to p = -l diffraction order which is centered at 
y = -MA ) and has two components (m, n) = (-1, 0) and (0, -1). The p = 0 
diffraction order, centered at y = 0, is a multiple-beam interference. 

All the different diffraction orders appear at the same time in the image 
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plane with different y-coordinates. The p - -1, 0, 1 diffraction orders 
are shown in Fig. 5.2. 

Two gratings of different period can be used. If a composite 
grating consisting of two gratings of small period difference is 
oriented along the y-axis, and each individual grating orientation makes, 
an angle 9 with the y-axis, but in opposite directions, the field of p = 1 
diffraction order, as in Eq. (5.7), is 

U(x,y; p = l) = A 0 B 1 U Q (-Mx sin 0 , -My - X d" 1 cos 0 ) 

+ A 1 B 0 U q {-Mx - X fj d” 1 sin 0 , -My - X f x d‘ l cos 0 ) , 

(5.8) 

where dj^ and d^ are the periods of gratings and respectively. 

The amount of shift in x and y coordinates is X f^ (d ^ + d^ ) sin 0 , 

-1 -1 

and Xf^tdj - d^ ) cos 9 , respectively. 

For the special case, 8 = 0, gratings Gj and G^ have the 
same orientation. The amount of shear is Xf^(dj^ - d^ ) directed 
along the y-axis; this is the smallest amount that can be achieved by 
two gratings of different periods. 

Two circular gratings of different frequencies having their 
center points in coincidence can also be used. The result is similar to 
the case of two linear gratings of different frequencies and has the 
same orientation discussed above. In this case, the shear is fixed for 
a pair of circular gratings. It is a circular shearing interferometer 


with constant radial shear. 



Figure 5.2. The field in the output image plane of the grating 
interferometer. Here we only show p - 1, 0, -1 
diffraction orders, and in p = 0 diffraction order 
we only show (-1, 1 ), (0, 0) and (1, -1) three diffrac- 
tion orders. There are only two diffraction orders 
in both p - 1 and -1 diffraction orders. 
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5.3. Variable, Purely Sinusoidal Wave Generator 

Purely sinusoidal amplitude gratings are very difficult to con- 
struct. By the use of the moire effect and spatial filtering, we can 
easily generate a purely sinusoidal wave. 

The optical setup is shown in Fig. 5. 3. Two superimposed 
gratings G 1 and G are placed at the object plane of a telecentric 
image-forming system. The field at the spatial-frequency plane is a 
collection of peaks with their locations distributed like the Miller's 
indices. A binary filter F is placed in the spatial-frequency plane. 

If the filter passes only two diffraction orders, e.g. (0, 1) and (1, 0), 
then the field at the output plane is 

U(x, y; z = 4f) = A^Bj exp(4-rrixd 1 sin 0 ) + A^q exp(-4TT i xd ^inQ) 

= | AqBjJ cos{4irxd 1 sin 0 +cp) , (5.9) 

if ] AqBj] = | AjBqI , or gratings Gj and G are identical, where cp 

-1 

is a constant. The sinusoidal wave has the spatial frequency 4 d sin 0 . 

Any two diffraction orders can be used for sinusoidal wave gen- 
eration. If the absolute amplitudes of those two diffraction orders 

I A B 1 and j A , B are equal, the generated sinusoidal wave 
! m n 1 1 m n 1 

2 

I(x, y) = | U(x, y; 4f)| has highest contrast. If the absolute amplitude 


sjt 2 

The area of each peak is roughly (X f) /A, where A is the area of 
grating being illuminated by the illuminating plane wave. 




Z = 0 2 f 


Figure 5. 3. Variable, purely sinusoidal wave generator. 
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strengths are not equal, the contrast will be reduced. We can vary the 
spatial frequency of the generated output sinusoidal wave either by vary- 
ing the angle 2 0 between gratings Gj and G^ , or by passing different 
diffraction orders in filtering plane. 

The output image plane is conjugate to the plane where the grat- 
ings are placed; hence the system is achromatic, although the spectral 
width of the light source is limited by separation of the adjacent diffrac- 
tion orders. A polychromatic light source can be used for illumination. 

5.4. Applications 

Most methods for optical information processing are based on 
interference. Interferometry is one of the most important methods in 
optical information processing. The grating interferometer derived in 
this chapter can be used in optical signal processing. New methods for 
optical testing and phase object visualization, image contrast reversal, 
spatial differentiation, complex amplitude addition and subtraction, 
image multiplexing, and fabrication of spatial filters will be presented 
in this section. 

5.4. 1. Optical Testing and Phase-Object Visualization 

The Ronchi interferometer has been used extensively in optical 
testing, since it is simpler and cheaper than most of the other methods. 
But the Ronchi interferometer has some drawbacks. To overcome the 

(4) 

drawbacks of the Ronchi interferometer is one of the motivations for 
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our new interferometer. Optical testing and phase object visualization 
can be achieved by the interference of the tested object wavefront with 
either a standard wavefront or a replica of the wavefront of the tested 
object with a shift. The grating interferometer can produce 
variable shear, either total shear or differential shear, and hence can 
be used in phase-object visualization. The new method is simpler, 
cheaper, and much more flexible than most other existing methods. 

The setup as shown in Fig. 5. 1 will be used in phase object 
visualization. For simplicity, we will assume that lenses and L 

have the same focal length f, and then the magnification factor M is 1. 

5.4. 1. 1. Phase-Object Visualization Using the Grating 
Interferometer to Produce Differential Shear 

The tested phase object U^(x, y) = exp[i0(x, y)] is placed in the 
object plane of Fig. 5,1. The two superimposed gratings and G^ 

have the same period d. The p= 1 diffraction order, centered at 
y = -A , from Eq. (5.7), is 

U(x, y; p = 1) = A () B 1 U Q (-x + A x ,-y - A y ) 

+ A 1 B 0 U 0 ( - X - A x' -V- A y> 

- A 0 Bi (exp[i« (-x + -y - A )] 

+ exp [i0 (-x - A x , -y-A )+ict]}, (5.10) 

if gratings Gj and G^ are identical, where a is a real constant 
depending on the relative shift of gratings G. and G . 

X 6 
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For an object with dimension in y direction 6 ^ Ay , the 

diffraction orders p are separated. Equation (5. 10) describes a two- 

beam differential shearing interferometer. The amount of shear is 

2 A = 2 X f d~* sin 0 . 
x 

The intensity of the p = 1 diffraction order is 
I(x, y) = | U(x, y; p = l)j 2 

= JAqBJ 2 cos 2 |[0(-x + A x , -y- A y ) -0(-x - A x , -y -A y )-a3 • 

(5.11) 


If A is smaller than the typical finer details of the phase object 
tested, we can expand the phase variation function 0 (x, y) in a Taylor 
series and retain only the first two terms, i. e. , 


0(-x + A x , -y-A y ) - 0(-x, -y- A y ) + A x 0 x (-x, -y - A y ) , ( 5 . 12 ) 

where 


0 Jx, y) 

X 


3 0(x, y) 
3x 


(5. 13) 


Equation (5. 12) can be rewritten as 

I(x, y) = | A q B 1 1 2 cos 2 |a x 0 x (-x, -y- A y ) - j } 

2 

= | A B | { 1 + cos (2 A 0 (-X, -y-A )-a)} . (5.14) 

*0 1' xx y 

The phase variation of the tested object becomes visible. The inter- 


ference fringes show the contour of constant derivation of the phase 
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variation of the phase object along the x-axis. If ct = TT > and A^0 x « 1 
Eq. (5. 14) can be written as 

I<x, y) = | AqBJ 2 (1 - cos (A x 0 x (-x, -y - A y )) 

* I 1 A 0 B l' 2 A x " y ' A y> • (5 - 15) - 

The output intensity is pi’Opor tional to the square of the derivative of 
the phase object. 

If the input phase object to be tested is a perfect lens, Uq(x, y) = 
2 Z 

exp[ia(x +y )], then the output image intensity from Eq. (5. 11) is 

I(x, y) = | A q B x | 2 [1 - cos (2a A x x) ] . (5.16) 

The interference fringes are straight lines parallel to y-axis. The 
parameter a can be measured from the period of the inteference 
fringes. If the input object is not a perfect lens, the interference 
fringes are not all parallel to y-axis. We can use this method to meas- 
ure the focal lengths of lenses and test the aberration of the lenses. 

5.4. 1. 2. Phase-Object Visualization Using the Grating Inter- 
ferometer to Produce Total Shear 

If the input object U Q (x, y) consists of two non -over lapping 

components U.(x - A , y) and U_ (x + A , y), 

1 X b X 


U n (x,y) = U,(x - A , y) + U,(x + A . y) . 


(5. 17) 


The objects U^(x, y) and U^(x, y) have dimensions § x and 6 in x 
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and y directions respectively, and 6 * A and 6 £ A . The p - +1 

' x x y y 

diffraction order in the image plane is 

V 1 (x,y) = A 0 B 1 U 0 (- X + A x , -y - A y ) + U Q (-x - A x> -y-A y ) 

= A q B 1 [U 1 (-x, -y-A y ) +U 2 (-x + E A x , -y-A y )] 

+ A 1 B q [U 1 (-x- 2A x , -y-A y ) +U 2 (-x, -y-A y )] 

= * A 0 B 1 U 1 ( ' X, ' y_A y ) + A 1 B 0 U 2 (_X ' " y_ A y } ^ 

+ A lB 0 U i(-x-2A x> -y-A y ) +A 0 B 1 U 2 (-x+ 2A x , -y-A y ) . 

(5. 18) 

The typical wavefield of V^x, y) is shown in Fig. 5.4. The wave- 
fronts of objects U 1 (x,y) and U 2 (x, y) overlap and are centered at 
(x = 0, y - - A ) as can be seen from Fig. 5.4. The intensity in the 
area where and U 2 overlap, when A^B^ = A^B^ exp(ia), is 

I(x, y) = | A q B x | 2 |U 1 (-x, -y - A ) + e ia U 2 (-x, -y - A y )| 2 . 

(5. 19) 

If Uj(x, y) = exp(i0 ^{x, y)) is the phase object under test, and 
U 2 (x, y) = exp(i0 2 (x, y)) is the standard with which the object U^(x, y) 
is to be compared, the intensity I{x, y) will show the deviation of 
0 (x, y) from 0 2 (x, y). If 0j(x,y) equals 0 2 (x,y), I(x, y) is just a 


uniform field. 




ypical wavefield 
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Usually the standard object U 2 {x, y) is an object that creates a 
uniform field, or a plane wave. The plane wave standard wavefront 
can be created by simply passing the light through the position of object 
U 2 (x, y). When (x, y) = 1, then 

I(x, y) = 1 A 0 B/ { 1 - cos (0 (-x, -y - A ))} , (5.20) 

if a = TT . For objects with small phase shift 0 « 1 , 

l(x,y) * I ] AgBjj 2 0 2 {-x, -y - A y ) . (5.2D 

The output image centered at (0, -A y ) is proportional to the square of 

the phase of the phase object, if the phase object has small phase shift. 

This is equivalent to what we get from the conventional phase-contrast 
(5) 

microscopy. 

Total shearing interferometry can be achieved by using only one 
grating as a beam splitter, i. e. , we can use the Ronchi interferometer 
to produce total shearing for phase object visualization. The tested 
object Uj(x, y) and the standard wavefront U 2 (x, y) are placed in the 
object plane with their centers separated by a distance \f /d in the 
direction of grating orientation. The result is the same as we have just 
discussed. 


5.4.2. Spatial Differentiation 

Optical spatial differentiation may be desirable in optical signal 
processing. In the performance of this operation, the principal problem 
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encountered in a coherent optical data processing system is the fabrica- 
tion of the spatial filter. A compound filter consisting of an amplitude 
filter and a phase filter to obtain the required transmittance function ju) 
has been reported. ^ One serious difficulty encountered in the fabrica- 
tion of this compound spatial filter is that the dynamic range of the 
photographic film is limited between zero and one. The alignment of 

the amplitude filter and the phase filter might be difficult. Computer - 

( 7 ) ( 8 ) 

generated holograms 1 and double exposure composite holograms 

have been used to bypass these difficulties, Birefringement elements 

such as Savart plates and Wollaston prisms have been used as beam 

splitters in a differential- shearing interferometer to perform the 

spatial differentiation. 

Differentiation can be approximated by a difference equation 
with a differential step e in the direction of differentiation, if the 
differential step e is smaller than the typical finer details of the 
object to be differentiated, 

£ U 0 <x,y) * i [U 0 (x+e, y) - U 0 (x,y)] , (5.22) 

if e ^ h . where h is the size of typical finer details of the object 
x x 

in x direction. Or equivalently, in the Fourier domain 

*[£ V*-^] 3 7 V V 'U> Hi v) • ,5 - 231 

if 

Uq( v , |i ) = 0 , for | v | £ h^ 
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In the region | v 1 £ h” 1 , sin ( j v) = y v , so that a sinusoidal grating 
can be approximately used as the filter for spatial differentiation. 

The essence of spatial differentiation is the creation of bipolar- 
ity (positivity and negativity) and the shift of the object wavefront. The 
grating interferometer can produce differential shear and bipolarity, 
hence can be used in spatial differentiation. Computer -generated 
holograms and double exposure composite holograms for spatial differ- 
entiation all can produce fixed differential step e . The new method 
uses the grating interferometer, which can produce variable differential 
step c to fit the various kinds of objects. 

The optical setup is as shown in Fig. 5. 1. The object Uq(x, y) 
is placed in the object plane. The field of p=l diffraction order, 
centered at y = , from Eq. {5/ 10) is 

U(x, y; p= 1) = A Q B 1 U Q (-x + A x , -y-A y ) 

+ A 1 B 0 U 0 < - X - A x > ■ y ' A y ) 

“ A 0 B 1 <2A *> ^ U 0 ( -*’ - y - A y> • <5 - 24 » 

if 2 A = 2 X f d * sin 9 is smaller than the typical finer details of the 
x 

object. If we rotate the composite grating (G^ + G^) 90°, we can have 
spatial differentiation along the y-axis, and if the orientation of the 
composite grating is 45° to the x-axis, we can evaluate 9 Uq/9x - 
9U 0 /9y . 
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5. 4. 3. Complex Addition and Subtraction 
In optics , two functions can be multiplied by using the 
Kirchhoff boundary condition. When two films are superimposed, the 
resulting transmittance is the product of their transmittances. Addi- 
tion of two images can be achieved by double exposure. But subtraction 
and division are difficult. Our grating interferometer can be used for 
real-time complex addition and subtraction. 

Many methods for subtraction using only optical elements have 
been already described. They are all based on the interference of the 
input signals. When they interfere constructively, we have the addition 
of the input signals; if they interfere destructively, we have the subtrac- 
tion of the input signals. The first method utilizes interferometry, 

( 9 ) 

either coherent holographic interferometry, total shearing inter- 
ferometry, or incoherent holographic interferometry. The 

second method uses a previously prepared hologram as a filter in a 

( 12 ) 

Fourier hologram recording geometry. The third method uses 

carrier -frequency photography to encode the polarity of each input sig- 
nal by shifting the modulation carrier -frequency half a period. 

The carrier -frequency modulation masks can be regular carrier- 
frequency masks or random diffusers. 

The methods that have been discussed either need very skillful 
laboratory technique or a 180° phase plate to achieve the bipolarity, or 
the grating must be shifted between exposures. Most of them are not 
real-time operations and are inflexible. 
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One method uses a prefabricated sinusoidal grating as a beam 
splitter in the Ronchi interferometer, ^ ^ The operation is in real time. 

The sinusoidal grating can be replaced by a Ronchi ruling. Our grating 
interferometer can be used as a total shear interferometer, since it can 
be used in addition and subtraction. By proper arrangement of the sig-_ 
nals to be added and subtracted, we can add and subtract four signals at 
once without using any phase plate, and the operation is in real time. 

For the addition and subtraction of two signals, one grating will 
be used as a beam splitter in the spatial-frequency plane of Fig. 5. 1. 
This method has been discussed in which a sinusoidal grating is used as 
a beam splitter. The input signals U^x, y) and U^fx, y) are 

centered at x = Xf/d and x = -Xf/d respectively, that is 

U 0 {x,y) = U x (x - \f/d, y) + U 2 (x +\f/d, y) . (5.25) 

The field centered at x = 0, y = 0 of the output plane is 

V(x,y) = AjUj (-x, -y) + A_ 1 U 2 (-x, -y) . (5.26) 

* 

The grating is real valued, hence A^ = A ^ „ Subtraction is possible 
if we shift the grating in the x direction, so that A = -A^ . Then 
Eq. (5. 26) can be written as 

V(x,y) = A 2 [Uj(-x, -y) - U 2 (-x, -y)] . (5.27) 

which is the subtraction of U 2 from . No phase plate has been 
used, and the operation of addition and subtraction is controlled by the 
shift of the grating, which is independent of the incident light 
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wavelength \ . 

For more than two signals, more gratings will have to be super- 
imposed for splitting the beams. For four input signals with dimensions 

in x and y coordinates smaller than min (A , A ), centered at 
7 x y 

(x = A , y = A ), (A , -A ), {-A , A ) and (-A , -A ), i. e. , the . 
x 7 y x y xy x y 

input object Uq(x, y) is 

U 0 <x,y) = U 2 (x +A x , y + A ) + U 2 (x + A x> y - A ) + U 3 (x - A x , y +A y ) 

+ U (x - A » y-A ) (5.28) 

4 x y * 

the output field centered at x = 0, y = 0 of the output image plane is 
V(x, y) = A Q B jU^-x, -y) +A 1 B Q U 2 {-x f -y) + -y) 


A o B iV 


-y) 


(5.29) 


The relative phases between A^, A B^, B ^ can be adjusted by 
shifting the gratings Gj and G 2 in their orientation directions. 

For adding or subtracting just two input signals, the use of two 
gratings as a beam splitter offers some advantage. From Eq. (5.29), 
if = 0, the output field is 


V (x, y ) = A n B 


A -1 B 0 

U l(-x,-y) +^-3- U 2 (-x,-y) 


(5.30) 


The distance between the centers of two input signals U. and XJ is 

A 

2 A x , which is adjustable. The ratio (A jBq)/(A q B can be adjusted 
by choosing different space-to-period ratio q of the binary gratings G^ 
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and , 


as can be seen from Eq. 


(4.2). 


5.4.4. Real-Time Image Polarity^Rever sal and Two- 
Dimensional Photometry- 

Several methods for complete or partial reversal of the polarity 

of images have been proposed. These are alternatives to the traditional 

photographic method. ^ ^ Image polarity can be reversed in carrier - 

frequency photography by using the nonlinearities of recording media, a 

combination of spatially modulated signals, or reconstruction from the 

zeroth diffraction order. It is also possible to reverse the polarity 

of an image by using carrier -frequency holography and making the signal- 

to-reference ratio greater than 1 in the recording of the hologram, or by 

(17) 

using the zeroth diffraction order in the reconstruction. All these 

alternative methods are two-step processes, and they are not superior 
in either convenience or image quality to conventional photographic 
polarity-reversal. To my knowledge, there are as yet no real-time 
optical methods for image polarity reversal. 

Image subtraction is the essence of image polarity-reversal. 

Our grating interferometer can be used in addition and subtraction, 
hence can be used for image polarity-reversal. 

In this method, a Ronchi interferometer is used as a total shear- 
ing interferometer. Only one grating is used as a beam splitter in 
Fig. 5. 1. Both the object Uq(x, y) and a uniform field C I(x - 1 A, y) 
of limited extent are placed in the input object plane. The centers of 
Uq(x, y) and the uniform field are separated by a distance t A = £ \f/d 
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in the x direction, normal to the grating slits, where ! is a nonzero 
integer. C is an adjustable amplitude level of the uniform field, 
j C | ^ 1, and I(x, y) is an indicator function, 

( 1 if ] x| £ 1/2 6 , and | y| £ 1/2 6 

I(x.y) = X Y 

( 0 otherwise , 

where 6 and 8 are the maximum dimensions of the object in x and 
x y 

y-directions respectively. 

The input object is U Q (x, y) + C I(x - iA, y). The total output 

field is 

U(x, y) = ^ [A n U Q (-x-nA, -y) +C A^ I(-x-nA , -y)] 
n 

= J V (x, y) , (5.31) 

n 

where 

V (x,y) = A U n (-x-nA, -y) + C A . I(-x-nA , -y) . (5.32) 

n n u n-Jt 

If the object size 6 is smaller than A, then for different n’s the 

x 

functions V (x, y) appear in the image plane separated in space and 

centered at (-nA, 0). The different diffraction orders V (x, y) have 

n 

different image -transfer characteristics. Here, the image-transfer 
characteristic is the relation between the transmittance of the input 
object Uq(x, y) and the absolute value of the output field amplitude 
j V (x, y)| . The image-transfer characteristic has the shape of a 
and intercepts the U^-axis at point b (if b > 0). The slope is -1 to 
the left of this point, and +1 to the right. The value b is the polarity 
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level. For output image lV (x,y)1 centered at (-nA, 0), from 

n 

Eq. (5.32), the polarity level b is -A^ (C ) *. From the defini- 

tion, it is clear that the image-transfer characteristic and hence the 
polarity level are functions of C and the A^'s. Therefore, the polar- 
ity level can be adjusted easily by rotating two superimposed polarizing 
sheets placed over either the object Uq(x, y) or the reference window 
I(x - f A , y). If the polarity level b is greater than max U^(x, y), the 
polarity of the object U^fx, y) is completely reversed, but if 
min Uq(x, y) < b < max U^(x, y), the polarity is only partially reversed. 
For 0 < b < min Uq(x, y), the polarity is not changed, but the contrast 
has been enhanced. 

Typical input and output for 1 - 2, = q • sinc(nq) • exp(inrt/2), 


as shown in Fig. 5. 5, are 

V (x, y) = A U (-x-n A , -y) + C A I(-x-nA,-y) , (5.33) 

Xi n u xx ** u 

Vq(x, y) = q t u 0 (-x, -y) - C sinc(2q) I(x, y) ] , (5.34) 

V x (x, y) = iq sinc(q) [U Q (-x- A , -y) - C I(-x- A , -y)j , (5. 35) 

V 2 (x, y) = q [-sinc(2q) U Q (-x-2 A , -y) + C I(-x - 2 A , -y) ] . (5. 36) 


The polarity levels for V Q , V^, and V are C sine (2q), C, and 
C/sinc(2q) respectively. The image-transfer characteristics for V^, 
Vj, and are shown in Fig. 5.5. We have outputs with different 
polarities at the same time. 



The image-transfer characteristics. 

(a) The polarity level b > 1, total image 
polarity reversal, (b) 0 < b < 1, partial 
polarity reversal, (c) 0 < b < min Uq , 
no polarity reversal, contrast enhance- 
ment. (d) b < 0, contrast reduction. 
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The input (the object U Q and the reference window) can be 
shifted through any desirable transverse displacement. The output 
image pattern will follow the movement of the inputs. 


5. 4. 5. Fabrication of Filters 

We have already shown that the grating interferometer can be 

modified into a variable, purely sinusoidal wave generator, and a 

purely sinusoidal grating can be used in spatial differentiation. In this 

section, we want to use the sinusoidal wave generator to fabricate 

spatial filters for some optical signal processing. 

Spatial filters for spatial differentiation, subtraction, and corre- 

, . (8,10,18,19) 

lation have been fabricated by the multiple -expo sure technique. 

2 2 

For example, in the fabrication of the Laplacian operator 9 /dx + 

2 2 

9 /9y , five exposures are made. A phase plate or movement of the 

recording film between each exposure is needed in order to produce the 
bipolarity. Film nonlinearity can create troublesome crosstalk, and 
the mechanical shift between each exposure is undesirable. 

All the filters fabricated by the multiple -exposure technique can 
be fabricated by the use of the sinusoidal wave generator. Only one 
exposure is needed. The relative sign of each point source in generat- 
ing the sinusoidal gratings can be controlled easily by the shift of the 
modulation gratings in their orientation direction. No phase plate is 


needed, and the mechanical movement between exposures is eliminated. 
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Several spatial filters for mathematical operations such as 
8/8x, 8 Z /0x^, 8 Z /8x^ + 8^/8y^, 8 Z /8x8y can be fabricated. We 
observe some mathematical identities that 

— = lim [iu(x, y) * f (x, y)l 
d h 0 

where 

f ^ (x, y ) = 6 (x + h) - 6 (x) , {5. 37 ) 


and * is convolution and 6 is the Dirac delta function; 


8u , 8u 
9x 8y 


lim 
h -» 0 



U (x, y) * f 2 (x, y) 



where 

f 2 (x,y) = 6 (x + h, y) + 6 (x, y + h) - 2 6 (x, y) ; 


2 2 

8 u , 8 u 

2 2 
8x 8y 


lim [4 U(x, y) * f (x, y) 1 , 

h -* 0 L h 3 J 


where 


(5. 38) 


f 3 (*. y) “ 6 (x + h, y) + 6 (x -h, y) +6 (x, y +h) + 6 (x, y - h) - 4 6 (x, y) ; 


8 


2 

u 


8x9y 


= lim _L U(x, y) * f (x, y) 
h 0 ,2 * 

XI 


(5. 39) 


where 


f 4 (x, y) = 6 (x +h, y +h) + 6 (x, y) - 6 (x +h, y) - 6 (x, y +h) . (5. 40) 


By the use of the convolution theory of Fourier transform, we know that 
the spatial filters for mathematical operations 8/8x, 9/8x +8/8y. 
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8 Z /9x 2 +9 2 /0y 2 , and a 2 /8xay are f j , f 2 , f 3 , and X^ respectively, 
and the filters are composed of sets of gratings. 

The spatial filters can be recorded holographically. The optical 
setup is shown in Fig. 5. 6. The left portion of the setup is identical to 
the sinusoidal wave generator in Fig. 5. 3. The field in the filter plane 
is a set of peaks (approximated to delta functions), 

u F (x,y) = A n B m 6(x " ^x’ y ‘ mA y } ’ (5 ' 4l) 

n m 

-1 - 1 

where A = Xf d sin 0 , and A = X f d cos 9 as defined before, 

x y 

We can adjust the angle 9 and the orientation direction of the composite 
grating to obtain the desirable array of peaks. The relative signs of 
different diffraction orders can be adjusted by shifting the grating along 
its own orientation direction, and the amplitude A^B ^ can be con- 
trolled by using proper space-to-period ratio q of the grating. We use 
a binary filter that passes only diffraction orders yielding either fj , 
f , f or f , and they will interfere with an offset plane wave to record 

M J 

a Fourier hologram. The recorded Fourier hologram is the spatial 
filter. 

The gradient correlation filter: 

If the input to the coherent data processing system is U(x, y), 
and the output of the system is 


g(x, y) = JJ grad U' (x, y) • grad U*(x'- x, y'- y) dx' dy' , 


(5.42) 
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which is equivalent to 

g(x, y) = lim I --r [6 (x + h, y) + 6 (x - h, y) + 6 (x, y +h) + 6 (x, y - A ) 
h+ 0 ( li 

-45(x, y)]* U 0 (x, y) * U Q (x, y) 

= lim I ~ f (x, y ) * U n (x,y) 1 , 

h 0 ( h 2 5 ° ) 

where 

fj.(x, y) = [6{x+h, y) + 6(x-h, y) + 6 (x, y +h) + 6 (x, y-h) - 4 6 (x, y)] 

* Uq(x, y) . (5.43) 

(Z 1 ) 

This gradient correlation filter is useful in signal discrimination. 

The subtraction correlation filter: 

If the input is Uj(x, y) and U^(x, y), the output is 
[U j (x, y) - U^(x, y)] * Uq(x, y). The required spatial filter is 

f^(x, y) = [6 (x, y) - 6 (x + s)] ^ U Q {x, y) , (5.44) 



where s is the separation between the centers of U^(x, y) and 


U 2 (x, y). 


Filters f^ and f^ can be fabricated holographically, similarly 

to the fabrication of filters f, , f_ , f_ , and f . . The diffraction orders 

1 Z 3 4 

pass filter F to interfere with the spectrum Uq(x, y) instead of an 
offset plane wave as shown in Fig. 5. 6. 
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5.4.6. Image Multiplexing 

In conventional carrier -frequency image multiplexing, as dis- 
cussed in section 4. 3, a grating is in contact with either the object or 
the film. The image of the object is an amplitude-modulated version of 
the object. Between each encoding of each image, the grating was 
rotated at angle 9 . Each object is encoded in different 0 -coordinates, 
the so called "theta modulation. " The modulation mask is usually a 
binary grating. The spectra of different objects occupy lines of different 
angles 0 , each split into many diffraction orders. Only one object can 
be recorded at that particular angle 0 . The number of objects that can 
be encoded on a single film frame is limited by the angle Cp occupied by 
the spectrum of each object. If we can use a perfect sinusoidal grating, 
then no high harmonics exist. The spectrum of each object will just 
occupy part of the space in any 0 -coordinate. And many objects can be 
encoded in the same 9 orientation; if we have sinusoidal gratings of 
different frequencies, then more images can be encoded into a single 
film frame. 

We have already shown that purely sinusoidal gratings of vari- 
able frequencies and 9 -orientations can be generated easily. This 
variable, purely sinusoidal grating generator can be used in image 
multiplexing . 

The optical setup is shown in Fig. 5. 7. A telecentric image- 
forming system images the grating plane into the object plane, and the 



O l +G z F OBJ IMG 



Figure 5. 7. An optical arrangement for variable-carrier-frequency photography. The 
frequency and orientation of the carrier frequency are variable by rotating 
grating Gj and G^ , and passing different diffraction orders in filter F. 
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object plane is imaged by a lens into the image plane. A binary filter is 
used in the spatial frequency plane to pass only two diffraction orders. 

A sinusoidal grating will appear in the object plane to modulate the 
object. The contrast and the frequency of the generated sinusoidal grat- 
ing can be controlled easily as described in section 5. 3. The system is 
achromatic and broadband source can be used. 

In principle, the number of objects that can be encoded into one 
single film frame is determined by the resolution of the grating and the 
size of the spatial spectrum of object encoded. We can use the generated 
sinusoidal gratings of different angles 0 and different periods to modu- 
late the encoded images in such a way that the modulated spectra of the 
images do not overlap and can be retrieved separately. In reality, film 
nonlinearity and limitations on the dynamic range of the film could be 


troublesome. 
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Chapter 6 


Variable Shearing Interferometry Based on Carrier-Frequency 
Photography - The Carrier-Frequency Interferometer 

6, 1. Introduction 

Carrier-frequency modulation and demodulation have been very _ 
successfully applied in many branches of electronic signal processing. 

In communications, the signals to be transmitted are real-valued. 

Many modulation schemes have been used to transmit real-valued func- 
tions, such as pulse -position modulation, pulse -amplitude modulation, 
and pulse -width modulation. Only one kind of modulation is needed to 
transmit a real- valued signal. In optics, the signals in general are 
complex -valued. A modulation scheme for complex signals takes twice 
the effort of that for real-valued ones. Diffraction gratings can be used 
as the modulators. For a pulse -modulated complex signal, the ampli- 
tude of each pulse is proportional to the amplitude, and the shift of the 
position of each pulse from the sampling points is related to phase of 
the signal at the sampled points. It takes both pulse -amplitude modula- 
tion and pulse position modulation to register a complex-valued function. 
This can be done optically by carrier -frequency holography or carrier - 
frequency photography. 

In carrier-frequency photography, a grating is imaged through 
the object, which is interposed between the grating and the film as 
shown in Fig. 6. 1. The image will register information about both the 
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Figure 6. 1. Different optical setups for double-exposure carrier 
frequency interferometer, (a) Using one lens (L^), 
(b) using two lenses (L^ and L.^), (c) using no lenses 
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amplitude and the derivative of the phase of the object along the orienta- 
tion of the grating. If the object is a pure phase object, the image is a 
distorted grating, and the distortion is related to the derivative of the 
phase of the object in the direction of grating orientation. The distor- 
tion of the grating is very difficult to detect. But if the distorted grating 
image is superimposed on a master grating, moire fringes will appear 
that reveal the information about the phase variation of the tested object. 
The encoding and decoding processes are shown in Fig. 6. 2. 

The effect of a distorted grating has been used by Rowland and 
Lord Rayleigh^ in their explanation of the grating ghost, that a lateral 
displacement of some grating grooves causes a deformation of a dif- 
fracted wave, i.e., a phase shift. The so called "detour phase" effect 

( 2 ) 

has been used by Leith and Upatnieks in holography and by Brown and 

(3) 

Lohmann in computer holograms. 

Carrier -frequency photography has been used in image subtrac- 
tion^ and phase object visualization, both by non-real-time 

double -exposure techniques. Like holographic interferometry, real- 
time signal processing can be performed if a reference wavefront is 
available. If the reference wavefront is the shifted version of the object 
wavefront, or a plane wave, the reference wavefront can be created by 
another grating which is identical to the grating used to encode the ob- 
ject wavefront. The second grating is used as a beam splitter and 


re combiner . 



Phase Variation 


Grating 


Image 


Encoding 



Image 



Decoding 


Figure 6. 2 


Encoding and decoding of phase object. 
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The encoding process can be done by an image -forming system 
or by the Talbot effect of regular gratings. Regular linear, checker- 
board, or circular gratings can be used in the Talbot interferometer. 

The Talbot interferometer as described by Lehmann and Silva' 1 ’ is 

for observing the derivative of the tested phase object, but not the phase, 
of the phase object itself, A random diffuser can be used as a modula- 
tion mask if an image-forming system is used. The new carrier- 
frequency interferometer can be used to observe either the phase or 
the gradient of the phase object tested. Natural white light and a larger 
source can be used when an image-forming system is used in the encod- 
ing process. 

The carrier-frequency interferometer can produce variable 
shear. It can be applied to visualization of a phase object, spatial dif- 
ferentiation, addition and subtraction, image polarity reversal and color 
complementarity, and equidensitometry. 

6. 2. Double -Exposure Carrier- Frequency Interferometry 

In the last section we have already shown that carrier-frequency 
photography can register both amplitude and phase information of ob- 
jects. In particular, phase objects can be visualized by carrier- 
frequency interferometry. In this section we will show the detailed 
mathematical analysis of double-exposure carrier-frequency interferom- 
etry. 
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Three different optical setups as shown in Fig. 6. 1 can be used 
in double-exposure carrier-frequency interferometry. A point source 
located at (x^y ') of the front focal plane of lens will produce a 

plane wave exp[-2TTi (xx / + yy')/Xf] to illuminate the carrier-frequency 
encoding mask (EM), The encoding mask can be a linear grating, a 
checkerboard grating, a circular grating, or even a random diffuser. 

The recording film is placed conjugate to the plane in which the encoding 
grating is located, i. e. , the grating will form an image on the film. 

The image of the encoding mask (grating) is created by the image form- 
ing lens as in Fig. 6. 1(a), a system of telecentric image-forming 

lenses and as in Fig. 6. lfb), or the Talbot effect of a regular 

grating as in Fig. 6. 1(c). All the different setups in Fig. 6. 1 have 
certain advantages and disadvantages. Spatial filtering is possible in 
the setup of Fig. 6. 1(b), but two lenses are needed. No lenses are 
needed when the Talbot-image property is used, but the requirements 
on the size and spectral purity of the illuminating source are more 
severe. 

The carrier-frequency interferometer is a shearing inter- 
ferometer. When linear gratings and checkerboard gratings are used 
as the encoding mask, it is a lateral shearing interferometer; when 
circular gratings are used, it is a constant- radial- shearing inter- 
ferometer. A random diffuser is used as the encoding mask because of 
the convenience when a large size encoding mask is needed. 
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6.2. 1. Linear Carrier -Frequency Interferometry 

The optical setup is shown in Fig. 6. 1(a). The encoding 

mask used is a linear grating with transmittance function 

G(x) = ^ A exp (2TT i nx/d). The film is placed in the plane conjugate 
n n 

to the plane in which we place the encoding mask, i. e. , d +d^ -f =0, 
where d^ is the distance between the grating and image-forming lens 
L_ . d. is the distance between the film and lens L„ , and f is focal 

2 l 2 

length of both lenses L^ and L^ . Two exposures will be made. In 
the first exposure, the object Uq(x, y) is placed at a distance D to 
the left of the grating. 

The field just before the grating is a tilted plane wave, 

U(x, y; -0) = exp[2ni (ax + P y) ] , (6.1) 

where 

a = -x'/Xf 

M-y'Af , 

where x' and y ' are the coordinates of the point source in the front 
focal plane of lens L^ . The reason for the oblique illumination will 
be clear later. The incident plane wave will be diffracted into many 
plane waves by the linear grating. The field just before the object 
Uj (x, y), as shown in Eq. B. 1, is 

U(x,y;D-0) = ^ C^exp |^2rri ] • exp(-2ni yy '/Xf) , (6.2) 

n 


where 
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C = 

n 


A n exp | 


r iTTD 

Uf 2 


, / 2 X /2 
(x +y ) 


J exp j^-i tt 


n 2 XP ~| 

,2 J 


exp 


2TTi n 


Dx' 
f 


]• 


(6.3) 


The field is made up of many tilted plane waves with direction cosines 

in x and y ^ ^ ^ ^ ) an d amplitudes C^ # Those tilted plane waves 

illuminate the object. The field Vjx/y; x',y ') in the image plane 

/ / 

created by the plane wave of direction cosines (| - yj , J f according 
to Eq. (A. 12), is 


V 

n 


(x, y; x ' y ') = A exp 
n 




2tt in 
Md 


(1 + M~ 1 )(x 2 + y 2 ) 1 • 


exp 


2n i xx ' + y r y~ 

. M Xf 



x 'P -y y 'p v 

f ’ M " f )> 


(6.4) 

if D is small enough that the integration in Eq. (A. 9) , can be per- 
formed by the method of stationary phase. M = d./d Q is the magnifica- 
tion factor. D has to be smaller than & 2 /X, where 6 is the size of 
the typical details of the object U 0 (x, y). 

The total field at the film plane is 


v ( x ,y; y') = J A n ex p[-j“ (1 +M _1 )(x 2 + y 2 ) 
n i 

exp(- 2rrinx/Md) U, (- ^ 

I'M d 


jexp 

"2ni 

(xx * 4 yy ' ) 

. M 

Xf 

x / D 

-y 

y r D \ 

f 

* M 

f / • 


] 


(6.5) 

Each diffraction order of the grating has produced an image on the film 
with lateral shift nXD/d. The intensity of the first exposure is 
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E 1 ( x , Y> x ', y') = |V(sc, y; x',y')[ 


-11 

n m 


* 

A A 
n m 


exp[- 2 tt i(n-m)/Md] ^ 


# / X 

* M-m 


m X D 


, nXD x / D -y 

+ _ 5 — r - m 

x / D -y y'D\ 
f * M " f /' 


( 6 . 6 ) 


In the second exposure, the second object U\ (x,y) is used to replace 

Ld 

Uj(x, y), and the grating is shifted laterally by a distance s in the 
direction of its orientation, i. e. , G^Xj ) = G(x + s ) =^T B^expUTTinx/d), 
and exp{ZTrins/d). The intensity of the second exposure is 


E^y.x'.y') B B exp[-2TT i(n-m)/Md] 

n m 

V TT / * 4 . n x /p -y y 'D \ TT * / x ^ m X D X 'D -y y 'D v 

2\ M d " f ’ M " f / U 2 \ M + d ~~T~ • M ~ ~T~ I * 

(6. 7) 

The total exposure is 

E(x,y; x\y ') = E 2 + E 2 = ^ I (x, y) exp(-2nipx/Md) , (6.8) 

P P 


•where 
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p <X ’ yl = I A m+p A m [ U 1 (' Si + ^ -9-1-^) 

m 


* / x 4 m X D x 'D 


* I x 

XU l(“M + d 


x D -y y D\ 
f * M " f > 


* exp[2rr ipe/d] U 2 (-* + ( 2^_ D - =£ , jj - ) 


X u 


* / x ( m X D x 

2 rM + “a - 


^ ~y ,y^D\] 

f * M f / J ’ 


(6.9) 


The only effect introduced by the oblique plane wave illumination 

is the shift of the image by x'D/f and y'D/f in the x and y directions 

respectively. The size of the illuminating point source, A , has to be 

s 

small enough that A g « f6/D. From now on, we will assume normal 
plane wave illumination, x ,= y'= 0, and Eq. (6.9) becomes 


m 


— + 

m XD 

11 

M 

d ’ 

M 

— + 

m XD 

-y 

M 

d » 

M 


+ exp(2TTipa/d)U,(-5 + ( ™.t.P> XD , 12 ] B * /. i + HlP IZ)] 

F * 2\ M d M ) 2\ M d * M /J* 

( 6 . 10 ) 

If the film is developed in the linear region of the T-E curve 
of the film, the transmittance of the developed film is 


T(x, y) = t Q + t x E(x, y) , 


( 6 . 11 ) 


where t Q and t } are constants, and E(x, y) = E(x, y; 0, 0). 
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To extract the information, the developed film is placed at the 
input object plane of a telecentric image-forming system. Each diffrac- 
tion order I will form an island in the spatial-frequency plane of the 
P 

telecentric image-forming system. If we use a binary filter to pass 
only one diffraction order, e. g. the pth, the field at the output plane is 

f(x, y) = -y) . (6.12) 

The output field is a function of the distance D between the grating and 
the object. The output field f(x, y) is the interference of object 
U (x, y) and U (x, y). The object U (x, y) can be a uniform field, 

1 L 6 

U (x, y) = 1, or a shifted version of U (x, y), U (x, y) = UL (x - e , y-e ). 
z x £ _L x y 

The carrier -frequency interferometer can be a two-beam or 
multiple-beam interferometer, depending on the A^'s. If the setup in 
Fig. 6. 1(b) is used for recording, spatial filtering is possible, and we 
can eliminate the unwanted diffraction order A 

m 

This interferometer can be applied to the visualization of phase 
objects, spatial differentiation, image polarity reversal, and subtrac- 
tion of images. 

6.2.2. Checkerboard Carrier -Frequency Interferometry 
Checkerboard (or cross) gratings can be used to replace the 
linear grating as encoding masks. In a linear carrier -frequency inter- 
ferometer, the direction of shear is fixed by the orientation of the 
grating ; the checkerboard grating can give more directions of shear. 
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A checkerboard grating can be represented as 

G(x, y) = / / A exp(2ir i nx/d__ ) exp(2ni my/d ) , (6. 13) 

' *—i urn ^ y 


n m 


where A is the Fourier coefficient, d„ and d„. are periods of the 
nm x y 

grating in x and y directions respectively. 

For simplicity, we assume this grating is illuminated by a 
normal incident plane wave. Two exposures are made, one with object 
Uj(x,y), the other with object U^(x, y). By the same calculation as in 
the last section, the total exposure of the film is 

E(x, y) = 12 C^^fx, y) exp(-2rr i mx/Md^) exp(-2-n i ny/Md ) . 


m n 


where 


n A A* lu l - + D H + < n+ P^ P \ 
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*V IT- 


H(-n 


+ (m+fj X D -y j_(n+p) X D 

M d • M d 

x v 


/ X 1 ID 

11 +P^ D \] 

\ M d 

M d / ( 

X 

y > 


(6. 14) 

where A ^ and are the displacements of the grating in x and y 

direction between two exposures. 

The film is developed in the linear region of the T-E curve. 



184 


and placed at the object plane of an telecentric image-forming system. 

A binary filter passes only (n, m)th diffraction order. The field in the 
image plane, from Eq. (6. 14), is C^^J-x, - y)' Shearing in both x 
and y directions can be obtained depending on the diffraction order 
(n, m) chosen in the filtering plane. Visualization of the gradient of 
phase objects in many directions and with different sensitivities is 
possible from one recorded film. 

6.2.3. Circular Carrier-Frequency Interferometry 

Linear gratings and checkerboard gratings can produce lateral 
shifts as shown in the last two sections. For circularly symmetric 
objects, radial shearing is more convenient for the interpretation of 
interference fringes. When circular gratings are used as the encoding 
masks, the interferometer is a constant-radial-shearing interferometer. 

Circular gratings are very similar to linear gratings, espe- 
cially at the region far from the center. A circular grating of period a 
can be written as 

. _ V /- 2tt i nr/a 

C (r ) " 2, C n 6 * < 6 * 15) 

n 

A normally incident plane wave illuminates the circular grating. The 
field at a distance Z behind the grating, as shown in Appendix A, 

Eq. {C. 8) is 

y 2 2 

U(r,Cp;Z)— / C exp(-irrXZn / a ) exp(2n i nr/a) , 
f f n 
n 


{ 6 . 16 ) 
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■which is very similar to the diffraction of a plane wave by a linear 
grating, as shown in Eq. (6.2). Using the same calculation as in the 
last two sections, we find that the field at the image plane is 

V(r, cp) = 2 C n exp [j^dT (1 +M_1 )r 2 jexp(2TTi nr /a) u 1 (fj + 

(6. 17) 

where is the object placed at a distance D right after the circular 

grating. 

The total exposure from double exposures, one with object 
U 1 (r,cp) and the other with object U 2 (r,cp) , is 


E<r,cp) 


-I 


C' (r , cp ) exp (2 rr i r p / Ma ) 
P 


(6. 18) 


where 


C'(r 


■ *> ■ I W»{ u ‘ ( s + ^ 


*/ r , mX D 


+ U (— + ( m+ P^ D _cp\ U*(— + 
2VM a ' v / 2\M 


, -cp 


(6. 19) 


We can put the developed film in the object plane of a coherent signal 
processing system, and pass only the diffraction order P. Then in 
the output image plane the field is (r, -cp ). The interferometer is 
a constant radial-shearing-interferometer. The amount, pXD/a , is 
proportional to D and is variable. It can be used for testing lenses. 
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6.2.4. Random Carrier -Frequency Interferometry 

Random diffusers have been used in holography to illuminate the 
object when recording holograms. In holography, coherence between 
object wave and reference wave is essential, but as it introduces strong 
coherent noise, it is undesirable. Holograms made of diffusely illumir 
nated objects can be scratched, spotted with dirt, and even broken into 
pieces without serious loss of information. The diffused-object holo- 

(9) 

gram has less coherent noise, ’ but requires larger bandwidth; other- 

( 10 ) 

wise the annoying speckle noise will be introduced. 

Random carrier -frequency photography has been used in image 
subtraction, phase object visualization, image multiplexing, and image 
contrast reversal. Double exposure has also been used. Theoreti- 

cally, everything random carrier-frequency photography can achieve, 
can also be achieved by regular carrier-frequency photography, and 
probably better. But sometimes for convenience, random carrier- 
frequency is more desirable. 

Random carrier -frequency interferometry is especially good 

for simple objects, e. g. for the measurement of the transverse move- 

(13) 

ment of a rough surface. 

6.2.5. Lensless Carrier -Frequency Interferometry 


We have shown in the appendix that regular gratings can create 
self images. In all the interferometry we have discussed in previous 
sections, lenses are used to image the grating into the recording film. 
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We can also use the self-imaging property of gratings in carrier- 
frequency interferometry. 

The recording setup for the lensless carrier -frequency inter- 
ferometry is shown in Fig. 6. 1(c). A linear grating is illuminated by 
a normally incident plane wave, a recording film is placed at a dis- 
tance Z q to the right of the grating, and the object is interposed 
between the grating and the film. The distance between the object and 
the film is D. 

The field at the plane just to the left of the object plane, by 
Eq. (B. 1), is 

U(x, ytZ -D-0) = X A exp[-iirX(Z -D) n^/d^ ] exp(2TT i nx/d) , (6.20) 

o * n o 

n 

By the use of the Kirchhoff boundary condition, the field just behind the 
object is 

U(x, y; Z q -D+0) - ^ A^ expf-iir X (Z q -D) n /d^] exp(2TTi nx/d) U (x, y) , 

n 

( 6 . 21 ) 

where U^x, y) is the transmittance of the first object. When the wave 
propagates a distance D and arrives at the film plane, the field is 
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U (x, y; Z Q ) =ffu(x 1 ,y 1 ; Z q -D+ 0) exp{^ [(Xj-x) 2 + (y r y) 2 ]} dx^ 

( 6 . 22 ) 

y 2 2 

= / A exp(-irrX Z n /d ) exp(2n i nx/d) 

*->n o 

n 

/u i (x 1 ,y 1 ) exp|i^[(x r x + ^) +(y 1 -y) 2 ]}dx 1 dy 1 , - 

(6. 23) 

2 

If the distance D is smaller than 5 /X , the method of stationary 
phase can be applied in the evaluation of Eq. (6.23), and 

U(x, y; Z )=^A exp(-in X Z q n 2 /d 2 ) exp(2rri nx/d) Uj (x - , y) . 

n 

(6.24) 


6 is the size of typical fine details of the object. The intensity of the 
first exposure is 

2 


E^x, y; Z q ) = | U(x, y; Z q )| 


(6. 25) 


v -1 V' * 2 2 2 

y y A A exp[-irr X Z (n -m )/d ] exp[2rri(n~m) x/d] 
n m 


/ nXD \ rr * / m XD \ 

M x ‘~ - y ) u i y ) • 


( 6 . 26 ) 


After the first exposure, the grating is laterally shifted a distance A 

in the direction of the grating orientation, and the object IL (x, y) re- 

£ 

places U^(x, y) for the second exposure. The total exposure is 


E(x 


, y; Z ) = ^ I (x, y) exp(2iri px/d) , (6.27) 


o *—> p 

P 
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where 

y*»y> s 2 A m+p A m ex P[ inX Z o p(p+2m)/d 2 ] 

+ exp»ni pA,/d> U 2 (x- , y ) U~(x-51|2, y )] . 

(6.28) 

The total exposure for the lensless carrier -frequency interferometry 

is very similar to the one when lenses are used, except that the mag- 

2 2 2 

nification factor is 1, and a phase factor exp[inX Z Q (n -m )/d ] . 

6.2.6. Discussion 

Three different optical arrangements can be used in the double- 
exposure carrier-frequency. An image forming lens is used in the 
first setup as shown in Fig. 6. 1(a). A telecentric image -forming sys- 
tem is used in the second optical setup as in Fig. 6. 1(b), and spatial 
filtering to eliminate the unwanted diffraction orders is possible. This 
system is very flexible. The lateral shift of the image due to the tilted 
plane wave illumination, from Eq. (6.9), is x'D/f . If the illumination 
source is not a point source, every point source will create its own 
image, and they are shifted with respect to each other, the resulting 
image will be blurred. For good results, the size of the source B 
must be less than (f§)/(MD), where § is the size of typical fine details 
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of the object. If D = 0, the images resulting from different point 
sources will coincide, and a spatially incoherent source can be used. 

In the third optical setup, shown in Fig. 6. 1(c), where the self- 
imaging property of the grating is used, the coherence requirement for 
the Talbot effect has to be met, besides the requirement for using an 
imaging forming system. 

6.3. Real-Time Carrier -Frequency Interferometry 

Carrier -frequency interferometry can be performed by a double- 
exposure technique as shown in section 6.2. Two objects with trans- 
mittance functions U^(x, y) and U^(x, y) can be registered in a single 
film frame by two consecutive exposures. The interference pattern 
between U^(x, y) and (x, y) can be obtained by reconstruction from 
the developed film as in double -exposure holographic interferometry. 
The double -exposure technique is suitable for the situation when the 
objects Uj(x, y) and U^fx, y) exist in time sequence. 

There are some objections to the double -exposure technique. 
Film nonlinearities and grain noise may cause some difficulties, and 
the exact control of all parameters in both exposures is difficult. We 
have always assumed that the distance D between the object and the 
modulation mask is small enough for the method of stationary phase to 
be applicable. All these difficulties and restrictions can be removed 


by real-time operation. 
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In the double -expo sure process, if the tested object is a pure 
phase object, the intensity distribution of the first exposure, with the 
object present, is a distorted version of the grating transmittance func- 
tion used as modulation mask. To visualize the phase information, we 
make another exposure on the film with the object removed. The moire 
fringes resulting from the sum of the distorted grating image and the 
undistorted grating image will reveal the phase information of the 
tested object. This procedure is similar to the one shown in Fig. 6.2. 
The moire effect of the product of distorted and undistorted grating 
images will produce the same moire fringes in the first diffraction 
order as the moire effect of the sum of those two grating images. This 
is shown in section 4.2. The two grating images can be added by the 
double -exposure technique, while the super imposition of two grating 
images (distorted and undistorted) can be implemented in real time and 
will be shown in this section. 

By real-tine operation we mean that there is no time consum- 
ing process, e.g. developing a film, involved in the whole process after 
the system has been set up. We put the tested object in and obtained the 
result at the same time. 

The real-time decoding can be achieved by the use of another 
grating at the conjugate plane of the encoding grating Gj as shown 

in Fig. 6. 3. A telecentric image-forming system added to the right of 
the decoding grating G will provide the necessary filtering. The 
observation plane will be the plane conjugate to the plane where the 
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Figure 6. 3. The carrier -frequency interferometer. 
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object is placed. 

6. 3. 1. A Quantitative Analysis 

The standard optical setup for the real-time carrier-frequency 
interferometer is schematically illustrated in Fig. 6. 3 For simplicity, 
we will assume that the focal lengths of all lenses used in Fig. 6. 8 are 
f . The system consists of two cascaded telecentric image-forming 
systems. The object planes of the first and second telecentric systems 
are conjugate. The gratings and are placed in the first and 

second object planes, respectively. The object is placed at a distance 
D to the right of the grating Gj and forms an image at Z = 8f + D. 

Two spatial filters F^ and F^ are used. 

A normally incident plane wave illuminates the grating G^ . 

The light field after the grating G^ consists of many plane waves, 
which illuminate the object. At the plane Z = 4f where the grating G^ 
is located, the wavefront is the sum of many laterally shifted object 
wavefronts with a defocus of distance D. Each shifted object wave- 
front propagates in a different direction. If a filter F is used to 
pass only finite diffraction orders, then only those orders passed by 
filter Fj will reach the plane Z = 4f. The grating G^ further dif- 
fracts the incident wavefront into more diffraction orders. Those 
wavefronts traveling with the same direction cosine can be extracted 
by a binary filter F^ located at Z = 6f , and the final image is formed 
at plane Z = 8f + D. As the field at the image plane is the sum of 
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many laterally (or, with circular gratings, radially) shifted replicas of 
the object wavefronts, this interferometer is a shearing interferometer. 


6.3.2. Theory 

For mathematical simplicity, we assume that a point source 
located at the origin of the front focal plane of lens is used to 

illuminate the whole optical system. A pair of linear gratings, checker- 
board gratings, circular gratings, can be used as modulation and de- 
modulation masks. In this section, a pair of linear gratings is used. 

The field just before the object, from Eq. (6.2), is 

2 2 

TT/ ni - V A -iTtXD n /d 2ffi nx/d 9cn 

U(x, y; D-0) = / A e e . (6.29) 

' n 
n 

By the Kirchhoff boundary condition, the field just behind the object 
U Q (x, y) is 

2 2 

TT/ n xn\ - TT / * V A -irrXD n /d 2ni nx/d ,, 

U(x, y; D +0) = U 0 (x, y) 2 _, A e e . (6.30) 

n 


The wave propagates through lens and arrives at the first filter 

plane F^ . The field at the plane just before the plane F^ , by the 
Fourier transforming property of lens , as shown in Appendix A, 
is 


U(x.y, Zf-0) = exp <x 2 + y 2 > I (f { -f . £ ) 


2 2 

X exp {-i-rrXD n /d ) 


(6.31 ) 
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The phase factor in Eq. (6. 31) is due to the fact that that object is at 
a distance D to the right of the front focal plane of lens . If the 
fundamental spatial frequency of gratings and * 1/d, is higher 

than the bandwidth of the object Uq(x, y), Uq ( ^ ~ d ' \f ) be 

separated into separate islands for different diffraction orders n . 

A filter with transmittance function Fj(x, y; n) is placed at the n-th 
island of the filter plane . The transmittance function F^(x, y; n) 
may be 0 (opaque block), 1 (open window) or a pinhole. The field after 
the filter F^ is 


U(x, y; 2f +0) = ^ exp ["- (x Z +y 2 ) 

n n L \f 


e„O n W,W n (i-H,X 


) . (6.32) 


where 


W 


n{rrhfi)= c o(o-t-ff) F i (x -^ n| • < 6 - 33 > 


The field at the plane just before the grating G^ , by the Fourier trans‘ 
forming property of lens » is 

U(x,y; 4f-0) = ^ A r e" 2n 1 nx/d f W^, y^ 

n J 


x exp {' 05 [ ( x i +x " or ) + (jr i +y)Z ]} dx i dy 1 ■ 


where 


(6.34) 
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W^(x, y) = JJ ^ n ^ x 2 ,y z^ eXf> E 2rri ^ xx 2 +yy 2 ^ dx 2 dy 2 * (6.35) 


The second grating has the transmittance function 

G (x) = / B exp(2rri mx/d). The field just after grating G is 
c m £ 

m 


U(x,y; 4f+0) = £ 

n m 


-2tt i(n-m) x/ d 

A h> e 
n m 


x // w n ,x - y)e 5 tp {- 1 xD[( x+x i" !! d 2 ) + (y + y 1 ) 2 ]K 1 dy 1 • 

(6. 36) 

From Eq. (6. 36), we know that if the wave propagates a distance D 

it will form an image of W (x, y). It is clear from the geometrical 

n 

optics that the plane Z = 4f + D is the image plane for Z = D. The 
field at Z = 4f + D is 


U(x, y; 4f +D) = 22 
n m 


A B W ( -x + 
n m n \ 


m XD 
d ' 


-y 


) 


X exp[-2rri(n-m) x/d] • exp£ 


-inXD(n-m) 1 

,2 J 


■ 2 


I p (x, y) e 


-2tTi px/d 


(6.37) 


where 

I 

P 


2 , ,2 


V x>y) = 2 A „ B „- P w »(- x +2 r -y) e ' inXD p /d • < 6 - 38) 

^ n 
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If the second filter is a binary filter and only passes the p-th 

diffraction order, then the field in the final output plane Z = 8f + D 
is 

V(x, y)= I p (-x, -y) (6.39) 

= exp(-irrXD p 2 /d 2 ) ^ A B W (x + XD, y ). (6.40) 

*—> n n-p n \ d / 

n 

The final output field is a sum of many laterally shifted signals 
W^x + XD, yj with amplitude strengths A^B^ . 

6. 3, 3. Lateral Shearing Interferometer with Variable Shear 
The real-time carrier -frequency interferometer can be used as 
a two -beam or multiple-beam shearing interferometer depending on 
how many diffraction orders that the filter Fj has passed, or on how 
many A^'s are not zero. 

If the filter F^ is a binary filter which passes only the 0-th 
and 1-st diffraction orders, and filter F_ is also a binary filter that 
passes only the diffraction order p- 1, the field at the image plane 
Z = 8f + D, from Eq. (6. 40), is 

% 

Vjfr.y) = A 0 B_ 1 U 0 (x + ^ , y) +A 1 B 0 U 0 (x,y) , (6.41) 

2 

where an unimportant phase factor exp(-irrXD/d ) has been omitted. 

This interferometer is a two-beam shearing interferometer. 
The amount of shear, XD/d, is easily varied in real time. If two 
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gratings and G^ are identical, |A^j = | | for all n, and 

A^Bq = A q B ^ exp(ia); the two beams have equal strength. a is a 

constant depending on the relative lateral shift in the x-direction. It 

can be adjusted so that a - TT , and then A^B^ = -A^B 

If the filter passes only diffraction orders 0, 1, and -1, 

and F passes only diffraction order p = 0, the final output field is 
6 


V 2 <*■ y> = A. jB _ ,U 0 ( x + , y ) + A 0 B 0 U 0 (x, y, + A[ B , V Q ( x - , y ) . 


XD 


(6.42) 


This is a three-beam interferometer. It can be used in high-order 
differentiation as will be shown in the later section of this thesis. 


6. 3. 4. Real-Time Image Holographic Interferometer 
If the filter F^ passes only the +1 diffraction order and has a 
pinhole at the center of the island of zeroth diffraction order, and F^ 
passes only the p-th diffraction order, the field in the output image 
plane is 

V(x, y) = aA Q B _^+ A i B 1 _ p U o( x + IT" » Y ) * (6* 43) 

where 

a = // u o* x> dxdy 

is the total transmittance spatially integrated over the whole object. 
For p = 1 

V(x,y) = aA 0 B_ 1 +A 1 B 0 U () (x,y) . 


(6. 44) 
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This interferometer is a two -beam total shearing interferometer; the 
output field is independent of parameters X and D . This optical sys- 
tem is achromatic, and natural white light can be used. 

The real-time image holographic interferometer can be applied 
to image contrast reversal and phase object visualization. 

6. 3. 5. Real-Time Carrier -Frequency Interferometry with 

Checkerboard Gratings, Circular Gratings, and Random 
Diffusers as Modulation and Demodulation Masks 

All the modulation masks used in the encoding process of double- 
exposure carrier -frequency interferometry can be used in real-time 
carrier -frequency interferometry. The optical setup is identical to the 
one shown in Fig. 6. 3 except that the gratings G^ and G^ are re- 
placed by checkerboard gratings, circular gratings, or random diffusers. 
The mathematical analysis and results are similar to those for linear 
gratings. Here we omit the details of the derivation and just present 
the results for checkerboard gratings and circular gratings. They have 
certain applications. 

For checkerboard gratings with period d in x direction and 

X 

d y in y direction the field at the plane Z = 4f + D is 

U{x, y ; 4f +D) = ^ zL S zL A mn B fp ex P[ -2 TTi(m-i )x/d x ~ 2rri (n-p)y/d ] 
m n i p ^ ' 

2 2 1 

i-rrXD(m-f) iTTXD(n-p) / , £ XD _lPXD\ 

2 - 2 J W mnl- X d x * “ y d y >’ 

x y 



(6.45) 



where 
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W 


mn 


(x, y) = If w mn^ X r y l^ exp ^ 2rri ^ XX l +yy l^ dx l dy l ’ (6*46) 


W mn* X ' y) = U Q (x, y) F J (x, y; m, n) , (6.47) 


and A and B„ are the Fourier coefficients of the transmittance 

nin i p 

functions of gratings G 1 and G , and Fj (x, y; m, n) is the filter 
transmittance function for the island of index (m, n) in the filter plane 
Fj . Equation (6.45) can be rewritten as 


U<x,y;4£+D) = J.'ZK (x,y)e- 2TTUx/d * e - 2 "Py /d y, (6.48) 
i P 


i,P 


where 




1^-irrXD 




A B W 

mn m-f,n-l mn 


(-x + ^ XD , -y+^\D), (6. 49) 

X y 

If the filter F^ passes only the (£,p) diffraction order, the output 
field is 


r <.p ( * ,y) = Vp ( ' x '' y ’ 


= e i««, P ) yy i x+ ^ XD y + ?^p, D ) 

^ mn m-i , n-p mn' d A » y d K > 
m n r x v 


(6.50) 


where cp(£,p) is a function of i and p only. 
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For circular gratings, if the filter passes only the p-th 

diffraction order, the final field is 

V(r,cp; p) = Y A B W (r - AD , cp ) exp(-iTTAD p 2 /a 2 ) , (6.51) 

Z— <r n n-p n \ d * 

n 

where a is the period of the circular gratings. 

The interferometer using checkerboard gratings can be used 
as a multiple -beam lateral shearing interferometer, and it becomes a 
constant radial shearing interferometer when the circular gratings are 
used. They can be used in spatial differentiation. 

6. 4. Applications 

The carrier -frequency interferometer is a variable shearing 
interferometer. It can be applied to phase-object visualization, image 
contrast reversal and color complementarity, spatial differentiation, 
and complex addition and subtraction. 

The spatial frequency of the gratings used in the carrier - 
frequency interferometer has to be higher than the bandwidth of the 
object spectrum in order to eliminate the aliasing and satisfy the 
sampling theorem. The size of the input object is not a limiting factor. 
But in the grating interferometer, the size of the object is limited to 

K 

Af/d, where d is the period of the grating used. This is the main 
difference between the grating interferometer and the carrier - 
frequency interferometer. 
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6.4.1. Phase-Object Visualization 

6.4. 1. 1 Phase-Object Visualization Using Double -Exposure 
Carrier-Frequency Interferometry 

Double -expo sure carrier-frequency photography has been used 
in phase object visualization. Linear gratings, ^ circular gratings, ^ 

(13) 

and even random diffusers have been successfully used as the 
modulation carriers. In this section we will only present the cases of 
circular gratings and linear gratings used as modulation masks to 
encode phase objects. 

The optical setup is as shown in Fig. 6. 1(a). For simplicity, 
we assume that M= 1, i. e. , d^ = d^ = 2f . Two exposures are made, 
one with the phase object under test U^(x, y) = exp{i0(x,y)) present, 
the other with the standard object U^(x, y) present. The standard 
object can be a plane wave, a known wavefront, or the shifted version 
of the phase object under test. 

When the modulation mask is a circular grating, the total 
exposure, as shown in Eq. (6. 18), is 

E(r, cp) = C ' (r, cp) exp(2rri rp/a) , 

P P 


where 



+ ( m+pUD 

a 



x m AD 
a 



+ u (r 4 . (ra-frlXD 
2 \ a 



+ m \D 
a 



(6.52) 
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The most often used standard wavefront is a plane wave, i. e. , 

U^x, y) = 1. Then before the second exposure, the tested phase object 
is just removed. Equation (6.52) becomes 


C' 

P 


(r, *> - ^ 
m 


C 

m+p m 



(m+p) \D 
a 


-cp 



m\D 

a 


. -cp 



( 6 . 53 ) 


Let us consider the situation when A - n\D/a is small, and 

the Taylor series expansion of the tested object U^(r - A , Cp ) = 

^i0(r A , cp) _ ex _ i r ^ > _cp) - A0 (r, -cp) + * • * j can be well approxi- 

L r J 

g 

mated by only the first two terms, where 0 (x* , cp ) = — 0{r,cp). Equa- 

r or 

tion (6.53) becomes 


C' 

P 


(r,cp) = 

m 


C C ( 1 
m m+p 


-ip A 0 (r , -cp) 
. r 

+ e 


(6. 54) 


ip A0 (r , -cp) 


= Cll + e 


= C e *p(^ P A0 r (r, -cp)) cos pA^(r, -cp)) (6. 55) 


where C is a constant. We can place the developed film in a tele- 
centric image-forming system and pass the p-th diffraction by using a 
binary filter in spatial frequency. The intensity of the output field at 
the image plane of the telecentric image-forming system is 


2 2 / 1 > 

I(r,cp ) = |C / (r,-cp)| = c cos p A^ir, cp)) 


= c [1 + cos (pA0 r (r, cp)) ] 


( 6 . 56 ) 
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The intensity reveals the radial derivative of the phase variation of the 
tested phase object. 

If the modulation mask is a linear grating, and the standard 

object is U (x, y) = 1, then the total exposure, from Eq. (6. 10), is 
c* 


E<*. y> = X I {x, y) exp(-2TTi px/d) , 


where 


y*.y> = X A m4p A m 

r m r 





+ mXD 



(6.57) 


We place the developed film in the object plane of a telecentric image- 
forming system to perform spatial filtering. If the system passes only 
the +1 diffraction order, the output field is 

V(x,y) = I p (-x, -y) 

= A o A i [ ex P( 2 TT is / d ) + u i( x + ^ • y ) u *( x > y)] 

+ [ ex P( 2 TTi s/d) + U^x, y) U*( x - ~ , yjj 

= AgA^exp^TTi s/d) + exp(i2 A 0 x (x, y)) ] , (6.58) 

if A = XD/d is smaller than the size of typical fine details of the phase 
variation, where 0 (x, y) = 9 0(x, y)/8x, and s is the shift of the grat- 
ing between exposures. The intensity of the output field is 
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I(x, y) = j V(x, y)| 2 = | AqAj ] 2 cos 2 (2 A0 x (x, y) + 2 tt s/d ) . (6. 59) 

The gradient of the phase variation, 0 (x, y), is revealed in the inter- 

X 

ference fringes. The result is similar to Eq. (5. 14) of the grating 
interf erome ter , 

Double -expo sure carrier-frequency interferometry can be per- 
formed via the use of an image forming system, or by using the Talbot 
effect. Lensless carrier -frequency interferometry can be used in 
phase object visualization; the result will be identical to Eq. (6, 56), 

and (6. 59) if the modulation mask is a circular grating or a linear grat- 

2 

ing, except that the contrast is reduced by a factor cos(XZTT/d ), 

where Z is the distance between the grating and the recording film. 

2 2 
If Z = 2 Nd /X , the Talbot distance, then cos(XZn/d ) = 1, and there 

is no loss of contrast; N is an integer. 

6.4. 1.2. Phase-Object Visualization Using Real-Time Carrier- 
Frequency Interferometry 

The real-time carrier -frequency interferometer can be used 
as a two-beam interferometer. For phase-object visualization, if the 
filter Fj is a binary filter which passes only the zeroth and first dif- 
fraction orders and F^ passes only first diffraction order, the field 
at the image plane Z = 8f + D, as shown in Eq. (6.41), is 

V(x, y) = A Q B_ 1 U Q (x + A, y) + A 1 B Q U 0 (x, y) . (6.60) 

We have assumed that Linear gratings are used. If gratings and 
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G. are identical, then I A | - lB | , for all n, and A n B - 

2 1 n' ' n' U -1 

AjBq e ia , a is a constant. The intensity of the output field is 
I(x, y) = | V(x, y)j 2 = 2 | AqB_ ^ | 2 £ 1 +cos [0 (x +A , y) - 0 (x, y)-a]3 

2 ] AqB j Z { 1 + cos [A0 (x, y) - a ]} (6* 61-) 

if A= XD/d is smaller than the size of fine details of the phase varia- 
tion 6 . The phase a can be adjusted by laterally shifting either 
grating or G^ . For a = Tt/2 , Eq. (6. 61 ) becomes 

I(x, y) = |A () B_ 1 | 2 sin 2 (y 0 x (x, y)) 

^ l A 0 B -l | 2 (6.62) 

if — A 0 (x, y) « 1. The output image is proportional to the deriva- 

c* X 

tive of the phase variation. The interference fringes can be inter- 
preted as the contours of the derivative of the phase object along 
x-axis. 

For observing the phase variation instead of the gradient of the 
phase variation of phase objects, the real-time image holographic 
interferometer can be used. If the filter F^ passes only the first 
diffraction order and has a pinhole at the center of the zeroth diffrac- 
tion order, and ¥ passes only the first diffraction order, the field 
in the output image plane, from Eq. (6.44), is 
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V a (x, y) = aA () B_. 1 + A 2 B 0 U 0 (x,y) , (6.63) 

The intensity is 

I a (x,y) = 1 V a (x, y)| 2 

= jA^^B jl sin 2 (-■ 0(x, y) ) f (6.64) 

if coefficients and B have relative phase difference fT , and 

a = 1 . It can be achieved if the images of grating and grating 

G are laterally shifted by half a period. The interference fringe is 

C* 

the contour of phase variation. As parameter \ does not appear in 
Eq. (6. 63), the system is achromatic; natural white light can be used. 
For an object with small phase variation 0 (x, y) « 1, Eq. (6.64) can 
be rewritten as 

l a (x,y) =“ I \ a 0 B _i y ) 1 2 * 

It is similar to phase -contrast microscopy. 

In the above discussions, we have always assumed that gratings 
Gj and G^ have the same orientation with their slits parallel to y- 
axis. For phase-object visualization using carrier -frequency inter- 
ferometry, the results are shown in Eqs. (6. 56), (6.59), (6, 62) and 
(6. 64). The interference fringes do not reveal the information about 
the sign of the phase variation (or gradient of the phase variation). 

The loss of sign results from the fact that what we observe is intensity, 
not the amplitude. 
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If a reference can be provided, the information about the sign 
of the phase variation (or gradient of the phase variation) can be pre- 
served. We rotate grating through an angle 9 with respect to 

y-axis, and rotate through an angle - 0. Following the same 

analysis as in section 6. 3, the field in the final output plane Z = 8f + D. 
of Fig. 6. 3, is 

2 2 ST 

V (x, y) = exp(-inXDp /d )/A B (x+(n-p)A ;y+(n-p)A ) 

p *—> n n-p x y 

exp[2rri (2n-p) y sin 9/d] , (6.65) 

where 

A = Xf cos 0/d and A = Xf sin 9/d 
x y 

If the tested phase object is placed in the input object plane of 
Fig. 6.3, and if the filter F^ passes only diffraction orders 0 and +1, 
and F^ passes only the +1 diffraction order, the field in the final out- 
put plane is 


V(x.y) = A ( ,B_ 1 U 0 (x-A x ,y-A > e 


- 2tT i y sin 9 /d 


, * D TT/ . 2rriysin0/d 

+ A i B o u o (x ’ y) 6 


( 6 . 66 ) 


and the intensity is 
I(x, y) = LV(x,y)| 2 



12 / 1 
cos ^tt y sin 0 /d + — A 


0 x (x,y)+- A # <*,y)) 


(6.67) 
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if A and A are smaller than the fine details of the object. In 
x y 

Eq. (6.67) a reference sinusoidal wave cos(ZTTy sin 0/d) serves as a 
reference. 

When the filter F passes only the +1 diffraction and has a 
pinhole at the center of the zeroth diffraction order, and passes . 

only the +1 diffraction order, the field in the image plane is 


V(x,y) = a o B 1 


-2TTiysin9/d, 
e 7 + 


A 1 B 0 U 0 


(x,y) 


2rri y sin 9 /d 


( 6 . 68 ) 


it is the interference between the tested object and a tilted plane wave. 
The intensity is 

i(x,y) = l a 0 B j \ 2 + l A 1 B 0 i Z+ 2 l A i B o a o B i I cos ( 4TT y sin 9/d+ 0( x ,y) + a ), 

(6.69) 

where a is a constant. For 9=0, Eq. (6. 69) is equal to Eq. (6. 64). 
6.4.1. 3. Experiment 

An experiment on double -exposure circular carrier- frequency 
interferometry has been performed. A helium-neon laser was used as 
a source whose beam was expanded to 50 mm. The circular grating 
used as the modulation mask had a period of 0.25 mm. The object was 
placed 1 mm behind the grating, and a camera was used to image the 
grating onto the film. An eyeglass lens was used as a test object. 

Figure 6.4 shows the result of testing the lens. Two bright fringes 
in the shape of concentric circles were observed. Another lens with a 
drop of plastic resin in its center was also tested. Figure 6. 5 shows 


the lens under test. 
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Figure 6.4. An eyeglass lens under test. The white fringes are 
fringes of equal radial phase derivative. The fine 
circular lines belong to the circular grating. 


Figure 6. 5. An eyeglass lens with a drop of resin in its center 
is the tested object. 
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6. 4. 2. Spatial Differentiation 

6.4.2. 1. Spatial-Differentiation Using the Double -Exposure 
Carrier -Frequency Interferometry 

Double-exposure carrier-frequency interferometry has been 
used in image subtraction, since it can produce bipolarity. Bipolarity 
is essential to subtraction. Differentiation is a special kind of sub- 
traction. 

The modulation mask (grating) can be placed either in the object 
plane or in the film plane of an image-forming system as shown in 
Fig. 4. 5. Incoherent light can be used and is preferred to coherent 
light for it eliminates the coherent noise. To achieve the bipolarity, 
the grating has to be shifted half a period between exposures. It is diffi- 
cult to move it exactly half a period. This mechanical movement of the 
grating can be eliminated if a polarization grating is used instead of an 
ordinary grating as the modulation mask. A polarization grating con- 
sists of mutually orthogonal stripes of polarizing material as shown in 
Fig. 6. 6. A polarization plate serving as an analyzer is placed in front 
of the grating. If the input light is unpolarized, natural light, and if the 
analyzer passes only polarized light in the x-direction, then the light in- 
tensity distribution in the film will be a binary grating consisting of 
black and white stripes. If we rotate the analyzer 90°, the grating 
image shifts half a period. By rotating the analyzer we can shift the 
grating half a period and create bipolarity. 

Two exposures are made. In the first exposure, we image the 
object U^^y) onto the film. In the second exposure, we rotate the 



I 


Figure 6.6. A polarized grating. Here f means polarized 
in y-direction, — *■ mean polarized in x- 
direction. 
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analyzer 90° and image the second object U^(x, y) onto the film. The 
total exposure is 

E{x, y) = I 1 (x, y) G(x) + I^fx, y) G(x - d/2) 

= exp{2TTi nx/d) [1^ (x, y) + e inTT I^(x, y ) ] , (6.70) 

where 1^ (x, y) and I^^y) are intensities of objects and . 

We can place the developed film in the input object plane of a coherent 
signal processing system and extract only the n= 1 diffraction order; 
the final image is proportional to 1 1 (x, y) - y)| , the difference 

between 1^ and . 

If the second object U^fx, y) is equal to (x - A, y), the 
first object with a shift A in the x-direction, as when the object 
Uj(x, y) or the camera set (film and grating) is shifted between expo- 
sures, Eq. (6. 70) becomes 

E (x, y) = ^ A n exp(2rri nx/d) [1^ (x, y) + e inTf (x - A , y) ] . (6. 71) 

The first diffraction order is proportional to 

Ij( x .y) - I(x - A , y ) — A —px* y) * 

if A is smaller than the typical fine detail of the object. 

Spatial differentiation along another direction can be obtained 
by shifting the object along that direction between exposures* Higher- 
order differentiation can be achieved by multiple exposures* The 
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polarity and intensity strength of each exposure is given in Eqs. (5. 38), 

2 2 2 2 

(5, 39), and (5.40) for the evaluation of 9I/9x + 9I/9y , 9 l/9x +9 l/9y , 
and 9^I/9x9y . 

6.4. 2. 2. Spatial Differentiation Using the Real-Time 
Carrier-Frequency Interferometry 

As the real-time carrier -frequency interferometer can 
produce differential shear, it can be used in spatial differentiation. 

The optical setup in Fig. 6. 3 is used. The object to be differ- 
entiated is placed in the object plane OBJ which is at a dis- 
tance D to the right of grating . The field in the final output plane 
Z = 8f + D, from Eq. (6.40), is 

V (x, y) = ^ A B W (x + \D y ) (6. 72) 

p * 7 n n n "P n ' ® > 

where 

W n( X+ if ’ y ) = U o(* + y ) 

if the n-th diffraction order passes the filter F^ , and 

W^|x + XD, y) = 0 , if the n-th diffraction order has been blocked 

by the filter F^ , 

For spatial differentiation, we can vary the filters F and F 

i c* 

for a particular kind of spatial differentiation. For example, we can 
have filter F pass only the diffraction orders n = 0 and 1, and F 
pass only diffraction order p = 1 . The output field is 
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Vj(x,y) 


= A 1 B 0 U () (x ) y) + A q B :i U(x + A, y) 

2 A i B o A el V x -y> • 


(6. 73) 


where A = XD/d . We have assumed A^B 1 = "^1^0 ’ aS can 
achieved by shifting the image of grating in the plane of grating 

G_ by half a period with respect to grating G . For example, if the 

L* Ld 

gratings G and G are Ronchi rulings, the transparent part of the 

-1 Ld 

grating G will image onto the opaque part of the grating G . There 

1 A 

is complete darkness after the plane where the grating G^ is located. 

For second order differentiation, the filter F^ passes only the 
diffraction orders n = -1, 0, 1, and F^ passes only the diffraction 
order p = 0. Then the field in the image plane is 

V 2 (x,y) = A_ 1 B_ 1 U Q (x + A, y) + A^iyx.y) + A^U^x- A, y) . 


A 1 B J A 


2 a 

ax' 


U Q (x,y) 


(6. 74) 


if 


A 0 B 0 


2 A i B r- 2 A.iB.i 


(6. 75) 


The equality in Eq. (6. 75) can be achieved by shifting either the gratings 
Gj and G^ and choosing gratings with the particular ratio q of space 
to the period of a binary grating that produces A^ = 4 2~ |A^ | . 
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For higher-order differentiation, we have to use checkboard 
gratings as modulation and demodulation masks. The field in the out 
put image plane, from Eq. (6.49), is 



*r)-H 

m n 


A B 

m, n m-f , n-p 


W (x + (p-m) A, y + (f-n) A) , 
mn 


(6. 76) 


if filter F passes only diffraction order (f,p) , where a multiplica- 

Lt 

tive constant has been omitted. W (x + (p-m) A, y + (X-n)A) 

mn 

equals Uq(x + (p-m) A, y + (f-n)A) if diffraction order (m,n) passes 
filter Fj , and equals zero if diffraction order (m,n) is blocked. 

For the evaluation of the Laplacian, the filter F^ passes dif- 
fraction orders (m, n) = (1,0), (0,1), (0,0), (-1,0) and (0,-1), and 
F passes only diffraction order (f,p) = (0, 0). The output field, from 
Eq. (6.76), is 


V(x,y) = A 10 B 10 U 0 fec - A, y) + - A) + y ) 

+ A -i,o E -i,o u o fe+ A - y ) + A o,-i B o,-i u o (x ' y + A) 


l 10 B 10 A ( 2 + 2 ^ B 0^ x ' 
\3x By / 


if 


3x 2 3y 2 


y) . 


(6. 77) 


A 10 B 10 = A 01 B 01 A -l, 0 B -1, 0 ~ A 0, -1 B 0, -1 ” “ 4 A 00 B 00 (6 * 78) 


The relation in Eq. (6. 78) can be achieved by shifting the gratings and 

choosing proper grating space- to-period ratio q. 

9 2 

For the evaluation of 0 t the filter F passes diffraction 

3x3y * 
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orders (m, n) = (1,1), (1,-1), (-1,1), (-1,-1), and F^ passes only 
diffraction order (i,p) = (0,0) . 

6.4.3. Real-Time Image Polarity-Reversal, Color-Complementation, 
and Two-Dimensional Photometry 

6.4.3. 1. Real-Time Image Polarity-Reversal and Two- 
Dimensional Photometry 

Image polarity reversal is equivalent to subtracting the 
image U Q (x,y) from a uniform background b to form 
b - Uq(x, y) . The constant b is the image polarity level. We have 
shown in section 5.4.4 that the grating interferometer can be used for 
real-time image polarity reversal and two-dimensional photometry. 
Here we will show that the real-time image holographic interferometer, 
discussed in section 6. 3.4, can be used for image polarity reversal 
and two-dimensional photometry. 

The optical setup is the same as shown in Fig. 6. 3 and is dis- 
cussed in section 6. 3.4. We will assume that the illuminating source 
is monochromatic and the object U^(x, y) is black, grey, and white. 
The field in the output image plane, from Eq, (6.44), is 

F(x,y) = C a Q (\) AqB j + A 1 B 0 U Q (x,y; \ ) , (6.79) 

where 

a Q (X) = fj u 0 ( x » y ; M dxdy (6.80) 

is the measure of the X spectral amplitude of the whole object (spati- 
ally integrated over the whole object). C is the transmittance of the 
pinhole. Uq(x, y; \ ) is the X -amplitude component of the light field at 
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point (x, y) of the object plane when the object is illuminated by a unit 
amplitude monochromatic plane wave of wavelength X . For a black- 
and-white object (no color), each X -component Uq(x, y; X) has the 
same amplitude distribution (as a function of x and y) as the object 
U Q (x,y) itself, i.e., if we disregard the color, Ugfc.y; X) is identi- 
cal to U„(x, y) for every X. The A 's and B ’s are the Fourier 
O’* n n 

coefficients of gratings and . Constants C, a^x), Aq and 

Bq are all real and positive (we have assumed the object Uq(x, y) is 
real and positive), and the phase of A^ and can be varied by 

shifting the gratings G^ and G^ along the direction normal to the 
grating lines. They can easily be adjusted so that A^ and B are 
out of phase; then Eq„ (6. 79) becomes 


V (x, y ) = A 1 B Q [b{ X ) - U Q (x,y; X) ] , 


where 


. (6.81) 


b(X) = 


CA o B -i 

A i B o 


ff u 0 (x »y » ) dxd y , 


(6.82) 


is the polarity level. For black-and-white objects, b(X) = B for all 
X , where B is a positive constant. 

It is clear that the output image will have the same color as the 
illumination monochromatic source. Positive images, negative images, 
or partial reversal of the grey scale of the object can be obtained, if 
we control the image polarity level. Image-transfer characteristics 
similar to those shown in Fig. 5.4 can be obtained. 
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The constant C can be controlled by placing a neutral -density 
wedge or two superimposed polarizers in front of the pinhole for real- 
time adjustment, C is confined in the region between 0 and 1. To 
increase the amplitude of the reference beam, we can increase a^ 
by passing some direct beam through Gj and U Q (x, y), and focusing 
it by lens into the center of the zeroth diffraction order in the 

filter F . 

6, 4. 3. 2, Real-Time Image Color- Complementarity and 
Polarity -Reversal 

In the last section we have shown that a monochromatic posi- 
tive image can be obtained directly from a black-and-white negative 
film. In this section, we will consider the situation when the input 
object is a color object and the illumination source is polychromatic. 

The real-time image holographic interferometer is achromatic 
and can be illuminated by a polychromatic point source. Here we will 
assume a ’white’ light spectrum W(X) uniformly distributed with 
density equal to a constant a between blue and red. Positive color 
images can be obtained directly from color negatives, if the illuminating 
source of the system is a white point source, and the color object itself 
or the modification of it is color-balanced. For a color-balanced ob- 
ject b(X) = B , a constant, i, e. , the object UQ(x,y) consists of the 
same amount of different colors. A black-and-white object is color- 


balanced. 



221 


For example, let a color object 'Uq(x, y) = Uq(x) be one- 
dimensional and consist of three equal regions of different colors. The 
first region is cyan (minus red), the second region is magenta (minus 
green), and the third is yellow (minus blue). The amplitude of the first 
two regions vary linearly, and the third region is constant as shown in 
Fig. 6. 7, where we have decomposed the object Uq(x, y) into three pri- 
mary colors, red, green, and blue. The amplitudes U^(x, y; X^), 

Uq(x, y; X^), and Uq(x, y; X^) are the amplitude components of the ob- 
ject light field when the object is illuminated by a unit monochromatic 
source of red, green, or blue light respectively. 

If the system is illuminated by a monochromatic light source, 
the object light field will be U^(x, y; X ), where X is the wavelength of 
the monochromatic point source. For example, if the monochromatic 
point source is red, the object light field is shown in Fig. 6. 7(a), if 
the input object is the object we have just discussed. Following the 
same analysis as in the last section, the intensity of the light field in 
the image plane, from Eq. (6.81), is 

I(x,y; = |V(x,y; X r >[ 2 


= IAqBJ 2 [b(X r ) - U 0 (x,y; \ r )f , (6.83) 


where 


b(X r ) = 


CA o B -i 

A 1 B 0 


II U o^' y; * r * dxdy 


(6.84) 
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b(X) is the X -component of the image polarity level, or the X- 
component of the reference beam. From Eq. (6.83) it is clear that 
polarity reversal for Uq(x, Y> ^ ) has been achieved. If bfX^) = , 

the output field |V(x,y; X^)| is shown in Fig. 6.8(a), The output image 
appears reddish. 

Because different X -components are incoherent with each other, 
when the object Uq(x, y) is illuminated by a polychromatic point source, 
we have to treat each X -component independently and then sum them up 
incoherently. Spatial coherence for each X- component is still main- 
tained, since the source is a point source. If the color object Uq(x, y) 

is illuminated by a polychromatic source of red and green light, then 

Z 2 

the image we observe will be |U 0 (x,y; X r )| + |U (x,y; X )| , which 

consists of three regions of green, red, and yellow. The color of the 
third region, yellow, is the incoherent sum of red and green. The image 
we observe is different to the original image. If this red and green 
source is used to illuminate the interferometer, the image we observe 
will be 

H*,y; X,) + Ife.y; X g ) - |V(x,y; X^ 2 + | V(x,y; X )| 2 

= |A 0 B . 1 | 2 ([ b a r )- u o ,x > y; X J-U be, y; X )] 2 } . (6.85) 

to y to 

In the example we have, as shown in Fig. 6. 8(a) and (b), we will observe 
an image consisting of three regions of reddish, greenish, and yellowish 
color. Comparing Fig. 6, 7 and Fig. 6. 8, we can see that the polarity 







Figure 6. 8. The three primary color transmittance of the output 
field )V(x)] when the reference strength is . 
Color-complementary and polarity-reversed image 
is obtained. 
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of each primary color has been reversed. That is, if at point (x, y) of 
the color object U Q (x, y) there is no green (or red), then at point (x, y) 
of the output image plane there will be a lot of green (or red), and vice 
versa (here we have assumed that the green- component of the reference 
beam is larger than max U Q (x,y; X g )) . 

If the object is illuminated by a polychromatic source of equal 
amounts of red, green, and blue, the image of the color object is the 
color object itself as shown in Fig. 6. 7. If we use this source to illumi- 
nate the interferometer, the output image intensity is 

I(x,y) = | V(x, y; X^)) 2 + |V(x,y; X g )| 2 + |V(x,y; X b )| 2 

= l A 0 B - 1 ! 2 t t b ^ r ) " u 0 ( x »y; X r )] 2 + l>(X g ) - U Q (x,y; x g) ] 2 

+ [t>(X b ) - U Q Cx, y; X b )] 2 } . . (6.86) 

The object we have, U^(x, y), has equal amounts of red, green, and 
blue; it is color-balanced, and b (X ) = B . For if b(X) = B = in 
Eq. (6.86), the output field |V(x, y; X)| is shown in Fig. 6.8. It con- 
sists of three regions of reddish, greenish and bluish color. The color 
of the output image is the complement of that of the input object U Q (x,y). 
A color- complemented and polarity-reversed image is thus obtained. 

For a polychromatic point source with white spectral density 
W(X), the X -component intensity of the image of the interferometer, 
from Eq. (6. 83), is 
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I{x ( y; X) = |V{x,y; X )( 2 = |A q B J 2 S(X) [B Q -U 0 (x,y; X )] 2 , 

if the object is color-balanced, and b(X) = B^ . For each wavelength 
X, the polarity of the point (x,y) at the image plane is reversed. 
Polarity-reversal in every X-component is the color -complementarity. 
If point (x,y) of the object passes only red light, then the image of that 
point from the interferometer will consist of all other colors except red 
(minus red), which is complementary to red. 

If the object is black-and-white, it is color balanced so the 
reference wave appears to be white. There is color dispersion in the 
first diffraction order. The period of the gratings and must 

be small enough to separate different diffraction orders after the color 
dispersion. The binary filter in the filter plane passes only 

the first diffraction order of all X wavelength components. Positive 
black-and-white images, or images with partial reversal of the grey 
scale of the input object can be obtained. 

If the object is red, even though the illuminating point source 
is white, the light after the object will appear red due to the spectral 
filtering effect of the object. The reference point source at the pinhole 
of the filter then appears to be red too. Image polarity reversal 

can be achieved, but not color-complementarity in this case. The final 
image is red but with contrast reversal. 

To obtain color -complementary images, the reference point 
source has to be white. In general, some adjustment in the input object 
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plane will have to be made to achieve this. This can be best understood 
by an example. Let us consider the object with its three primary color 
components as shown in Fig. 6. 7 but with the third (yellow) region elimi- 
nated. This object is not color balanced; blue light is dominant, so the 
reference point source appears to be bluish. To obtain a white refer- 
ence source, we have to modify the object by adding a third region of 
yellow color (minus blue) to achieve color balance. After we add the 
third region of yellow color, this modified object is then identical to the 
object shown in Fig. 6. 7. Following the same analysis as given before, 
a color- complemented and polarity- reversed image is thus obtained for 
the modified object and hence for the input object. Only optical ele- 
ments, such as gratings and lenses, are used. No complicated and 
expensive photochemical process or electronics TV technique is 
required. 

The contrast level can be changed independently of the color by 
varying the strength of the reference beam, b(\ ) = , as explained in 

the last section. To understand this possibility, let us consider a 
second example, in which we have adjusted the strength of the refer- 
ence beam so as to equal C as shown in Fig, 6. 7. The resulting 

£ 

image is shown in Fig. 6. 9. In this case, the image of the object is 
color complemented, but only partially polarity reversed. 







229 


6.4. 3.3. Comparison of Various Methods for Image Polarity- 
Reversal 

Two methods for real-time image polarity reversal have been 
discussed in this thesis, one using the grating interferometer, the other 
using the real-time carrier-frequency interferometer. In the first 
method, monochromatic light is necessary, but white light can be used 
in the second method. Spatially coherent light is necessary. The size 
of the object should be smaller than X f/d in the first method, and the 

dimension of the spatial frequency of the spatial bandwidth of the 

. -1 

object should be smaller than the frequency of the grating used, d , 
in the second method. 

The achromatic property of the second method enables us to 
obtain a positive black-and-white image directly from a black-and-white 
negative film, and a positive color image directly from a color nega- 
tive. Gratings of higher frequency are needed to ensure that the spatial 
spectra of different diffraction orders do not overlap after the color 
dispersion effect. 

(16 17 18 ) 

Previous methods 9 9 using two- stage carrier-frequency 

photography require a different carrier -frequency photograph to be 
prepared for each different output image polarity. Good skill in photo- 
graphic process control, or precise control of grating movement is 
necessary in those methods. A carrier-frequency photographic negative 
can be used as the input object in the real-time carrier-frequency inter- 
ferometer, and it will still provide images with continuously variable 
polarity from the same carrier -frequency photograph. 
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Since the polarity level can be easily varied continuously in real 
time, both methods discussed in this thesis can be used for two- 
dimensional photometry. ^ ^ 

This real-time image-polarity reversal technique is applicable 
to any image processing such as photography and xerography, photome- 
try, spectral analysis, interferogram analysis, and image evaluation 
in aerial photography, electron microscopy, and medical diagnosis. 

6.5. Focusing and Alignment 

In the previous sections, we have often had to place a point 
source in the front focal plane of a lens, and an object, a grating, or 
a random diffuser in the back focal plane. It is difficult to place them 
exactly at the position we want. New methods for focusing and align- 
ment will be presented in this section. 

The self-imaging property of a regular grating, or Talbot imag- 
ing, is a well known phenomenon. A point source is used via a lens to 
illuminate the grating as shown in Fig. 6. 10. If the grating is illumi- 
nated by a plane wave, the Talbot image is identical to the grating itself. 

( 22 ) 

Cowley and Moodie have shown that if the wavefront incident on the 
grating has curvature 1 /Z , the position of the self-image of the 
grating Z = Z^ is given by 

Z' 1 + Z" 1 = (NZ t ) _1 , (6.87) 

2 

where N is an integer, 2 = 2d / X the Talbot distance. The 
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Figure 6. 10. A Talbot image forming system for focusing 
and alignment. 
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image will be magnified by a factor 

M = 1 + Z /Z r . (6.88) 

The relation between the curvature of the incident wavefront and the 
magnification of the self-image has been used to monitor the beam 
collimation. By beam collimation, we mean placing a point source in 
the front focal plane of a lens, so that the wavefield after the lens is a 
plane wave, i.e. , Z ^ = 0 . To monitor the beam collimation we place 

c 

two identical linear gratings to the right of the lens, one after the other. 
If the incident wavefront to the first grating is not a plane wave, the 
Talbot image of first grating will have a period different from that of 
the grating itself by a factor , defined by Eq. (6.88). This self- 

image magnification can be detected by the moire effect of two gratings 
of different period. This method of beam collimation and its sensitivity 
have been analyzed. 

We have shown in Appendix C how the Talbot effect occurs for 
circular gratings. As Eqs. (6.87 ) and (6. 88) also apply to circular 
gratings, circular gratings can be used in beam collimation. A circular 
grating can register the position of its center. Owing to this property, 
circular gratings can be used in alignment. Two circular gratings are 
placed to the right of a lens as in Fig. 6.10. The coincidence of the 
centers of the second grating and of the image created by the first 
grating will show the alignment of two centers of the gratings and the 
point source. The collimation can be checked by the moire fringes 
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created by the second grating and the image of the first image. The 
sensitivity of this method is similar to what we can get when linear 
gratings are used. 

In beam collimation, we test the flatness of the incident wave- 
front. We want to place a point source in the front focal plane of a lens 
to create a plane wave. In certain situations, we want to place some 
objects, such as gratings, diffusers or spatial filters, in the back focal 
plane of a lens, i. e. , we want to put an object where the wavefront 
incident to it has infinite curvature, or Z^, = 0. If the grating is suffi- 
ciently close to the focal point, the area of the grating which is illumi- 

2 

nated by the incident beam is about pAz/f) , where f is the focal 

length of the length, B the size of the incident plane wave, and Az 

is the distance between the focal plane and the plane in which the grat- 

2 

ing is placed. If the area (B Az/f) is small, the wavefield illumina- 
tion of the grating no longer can be approximated by a infinitely wide 
wavefield. The Talbot effect has somewhat broken down. Equation 
(6.88) can be explained from the point of view of geometrical optics 
and will remain valid even when the Talbot effect has broken down. 

Experimentally we found that if the finer detail of the object is 
smaller than the size of the Airy disc Xf/B, which is the smallest 
spot size we can focus our incident plane wave into, we can always 
place the object in the back focal plane. We can identify the focal 
plane by shifting the object around the focal plane. We can see a 
speckle pattern at a distance Z from the focal plane, and the closer 
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the object is placed to the focal plane, the bigger the speckles are. 
The change of the size of the speckles is approximately inversely 
proportional to AZ . 
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Appendix A. Imaging and Fourier Transforming Properties 

of Lenses 

A. I General Input, Output Relationship 

The object Uq(x, y) is located at a distance Zq in front of a 
positive thin lens of focal length f. The tilted plane of illumina- 

tion is created by a point source located at coordinate (x*, y # ) of the 
front focal plane of lens L^ , as shown in Fig. A. 1. 

For simplicity, we will use only one -dimensional notation. It is 
straightforward to extend this discussion to the two-dimensional case. 
Coordinates x , x , and x refer to the object plane, lens plane, 

X Cj 

and observing plane respectively. The illuminating tilted plane wave 
can be expressed as 

U(x 1> y^ z = -0) = exp[ 2-rri (ax +0y)] , 

where 

a = -x'/\f , P = -y'/Xf (A. 1) 

By the successive use of the Kirchhoff boundary condition, and 
the Fresnel diffraction formula, the field at the observation plane is 

u i (x) = //w exp iax i ) exp [r^ <x r x 2 )Z ] exp [-Tr x 2 ] 

x exp [rf (x 2‘ x,z ] dx i dx 2 ■ < A - 2 > 
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where a multiplying constant is omitted,’ The consecutive phase factors 
in Eq. (A. 2) are due to the tilted plane wave illumination, propaga- 
tion from the object plane to the lens plane, the quadratic phase factor 
of the positive thin lens, and propagation from the lens plane to the 
observation plane. 

If we assume that the lens aperture is large enough, the integra- 
tion limits over variable can be well approximated by , 03 ) « 

Equation (A, 2) can be rewritten as 


U.(x) = exp(iTT x 2 /X z.) Jj U 0 (xj) exp ^2tt i aXj •+ x^ j 

X e*p[*2 +2*2 (~ dx l 


where 


(A. 3) 


- = — + — - -f 

a z. I 

0 l 


We can perform the integration over x by completing the square of 
the quadratic phase factor of x^ . By the Fresnel integration formula. 


we have 
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/ 2 . 

/ . iTTXj \ 

U q ( x ^ ) explZniax 1 + - J 

0 


X exp 




a Z b \ 2 1 | 2tt i ab . , , 

exp I ax - iTT b X 

L Vi 


'-[fcN- £>)■*] 

/ I W" P Fb *1 +b ( Xa *^;)] d.. 


where 


-1 -1 -1 ,-l / „ 

a - z + z. -f t^O 
0 i 


. -1 _ -1 -2 
b - z Q - a z Q 


It can be shown, after some algebra, that 


■1 


(A. 4) 


(A. 5) 


(A. 6) 


1 + 1 


z i z o - b f 


(A. 7) 


i. e. , if the object moves a distance b closer to the lens, then it will 
form an image at the observation plane. Or, we can say that the object 
is out of focus by a distance b. 

Equation ( A. 4) can be extended into two-dimensional cases 


easily. 
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U.(x, y) - exp | 1 
i 


a a b \ , 2 2. 

T 2 ' x +y } 

1 Z i Z 0 


exp /- TT — - b \(a x + p y ) - in b X (a. 2 + P ) 


\ z O z i 

/ u 0 ( x r y i ) j e xp p b | 


X l +b [ Xa z 




(A. 8) 


If we place our object at z^ = - D, and the observation plane 

is the image plane of the plane z^ - d^ , i. e. , z^ = , and 

d^ + d * = f 1 , then a = d^fd^ - D)/D , and b = -D . The field at 

z. - d. is 
l x 


j u o (x i*yi )ex p 



x-UD 
i 



■ (A. 9) 


i. e. , if the object field U (x^y^) propagates a distance -D from 
where it is located, then it will form an image at z. = d. # 

There are two interesting special cases: (1) a * = 0 , the 

observation plane is the image plane, and (2) = f , the observation 

plane is the back focal plane of the lens. 


A. 2 Imaging Property of Lenses 

If the observation happens to be at the image plane, then a ^ =0. 
We have no quadratic phase factor in the integration over variable x_ in 
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Eq. (A. 1.3). Integration over x 

M 

Subsequently an integration over x 
tion 


1 X 

yields a delta function 6 1 j — 

1 d o d i 


1 


yields, in two-dimensional nota- 


U.(x, y) = expj^-— - (1 +M 1 )(x 2 + y 2 )Jexp[-2TT i(ax + (3 y)/M]U 0 (- “ - ) , 


where 


(A. 11) 


M = d. /d . 

l 0 


The field at the image plane has been reversed and magnified M times. 
The introduction of tilted plane wave illumination yields a phase factor 
exp[-2rr i(d x y)/M] in the image field. When the image is detected 
by a square-law detector, such as film, the phase factor is not impor- 
tant. The angle of the illuminating plane wave has no effect on the final 
image. Consequently, point-source illumination is not necessary, and 
spatially incoherent illumination can be used. 

If the distance D in Eq. (A. 9) is so small that \D is much 
smaller than the typical finer details of object U^(x,y), then we can 
apply the method of stationary phase to perform the integration over 
Xj and y in Eq. (A. 9), and 

V(x, y) - exp[g- {1 +M“ 1 (x 2 +y 2 )jexp [~2ni{ax + 0 y)/M +iTTXD(a 2 +p 2 )] 


x u o(-^i +XaD ' * 


(A. 12 ) 
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The image has been shifted an amount X&DM in the x-direction and an 
amount X^DM in the y-direction, because of tilted plane wave illumi- 
nation. 

A. 3 Fourier Transforming Property of Lenses 

The other case of special interest is when z. = f - the observa- 
tion plane is the back focal plane of the lens. 

When z. = f, then a = z n , and Eq. (A. 1. 8) can be rewritten as 

l 0 ’ 

u i <x) = exp [rf (* -T) (xZ+y2) K(fj fr p )’ (A - 13) 

where 

U 0 (x.y) =/ U 0 (x^, y^) exp[2rri(x^x + y^y)] dx^yj^ . (A. 14) 

The field at point (x, y) of the back focal plane of the lerts is 
related to the Fourier spectrum of the object field at frequency 
(x/X f - a , y/Xf - (3 ). There is a quadratic phase factor, which will 
vanish if z^ = f , i. e. , if the object is located at the front focal plane. 
The tilted plane illumination introduces a shift of amount (aXf , $Xf) 
in the observation plane. Spatially coherent illumination is needed for 


spatial filtering. 



Appendix B. Spatial and Temporal Coherence Require- 
ments of Talbot Effect 


B. 1 Talbot Effect 

A monochromatic point source located at point (x / , y / ) of the 
front focal plane of a lens is used to illuminate the grating located to the 
right of the lens as shown in Fig. B. 1. The grating is illuminated by a 
tilted plane wave exp[-2rri (x'x + y / y)/Xf] which is created by the point 
source. The incident tilted plane wave will be diffracted by the grating 
into many plane waves, each of which propagates with different phase 
angle and phase attenuation. By the use of the Fresnel diffraction 
formula, the field at distance z beyong the grating is 

U(x, y; x ,y ; z; X) = ^A^ exp[-irT z(x'^ +y'^)/X f^ 


exp [-2rri(xx' + yy')/X f ] exp [-in X zn /d ] exp(2Trin(x/d +zx' /fd)] 


(B. 1) 

where G{x) = / A exp(2mnx/d) is the grating transmittance function, 
n 

X the wavelength of the monochromatic source. The first two phase 
factors are from the tilted incident plane wave, the third phase factor 
is the phase attenuation of the n-th diffraction order, and the fourth 
phase factor is the n-th diffraction order with a phase shift zx'/f due 
to tilted plane wave illumination. 

For simplicity, we assume that the grating is a pure sinusoidal 
grating, i. e. , A^ - 1, Aj - A^j = ~ . We are only interested in the 
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first harmonics of the image at distance z beyond the grating. The 
image is 

3 2 

l(x, y; x', y ' ; X ; z) = — +2 cos(n Xz v ) cos[2tt v(x + zx' /f)] 

+ — cos [4rr v (x + zx' /f) ] , (B. 2) 

where v = 1/d. The first harmonic is 

2 

Ij(x, y; x', y'; X ; z) = 2 cos ( tt X z v ) cos [2ttv <x + zx' /f) ] . (B. 3) 

2 

The constant factor cos(ttX zv ) is the reduction factor of the con- 
trast. Three cases are of interest: 

/ ^ Mv 2 /x = z T 

z= | (2M+l)v Z /X = Z c (B. 4) 

V (2M+2)v Z /X = Z q • 

They correspond to cos(Z^) = 1, -1, and 0. The image at Z = Z,^ 
is a replica of the grating itself, the image at Z = Z has contrast 
reversal, and the one at Z = Z^ has zero contrast. The image has 
been shifted by Zx'/f. This interesting feature of sinusoidal gratings 
holds also for Ronchi rulings. 

B.2 Spatial Coherence Requirement - Source Size 

If the source is not a single point, then each point of the source 
will produce a fringe. The fringes produced by each point source are 
shifted relative to each other, as can be seen from Eq. (B. 3). The 
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contrast of the resulting fringes will be" reduced if the source is too 
big. For a source of size B, the image in the first harmonic is 

B/2 

I(x, y;X;z) = / I (x, y; x', y' ; X ; z) dx y 

-B/2 

2 

= 2 B sinc(B vZ/f)cos(rrXZv )cos(2tt vx) . (B. 5) 

The constant sinc(Bv Z/f) is the reduction of the fringe contrast due 
to the broadness of the source. For good fringe contrast, 
sine (Bv Z/f) £ 2/tt , the size of the source B should be smaller than 
f /2 v Z , 

B £ f/2 v Z . (B. 6) 

B. 3 Temporal Coherence Requirement - Monochromaticity 

For different wavelengths the Talbot distance varies according 
to Eq. (B. 4). If the point source is not monochromatic, then color 
dispersion will appear. For good contrast, the spectral width of the 
polychromatic light must be limited. 

For a polychromatic source of spectral width AX , and mean 
wavelength X . the resulting image in the first harmonic is 

2 

I 1 (x,y;x / =0, y'= 0; z) = 2 cos(tt X V Z Z) sinc^^-^-| — — )cos(2ttvx) . 


(B.7) 
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Let Z 


= 2 


Md Z /X 


the Talbot distance for wavelength X > then 


1^-2 sine (MA X /X ) cos(2ttvx ). 


(B.8) 


For good contrast, MAX /X ^ 1/6 


(B. 9) 



Appendix C. Talbot Effect of Circular Gratings 


A circular grating C(r) = ^ C exp(2ni nr /a) is illuminated by 

n 

a normal incident plane wave. Using the Fresnel diffraction theory, 
we find the field at a distance Z behind the circular grating is 


U (r 


2tt 00 

,cp ; Z) = y* C / / exp{2rr inr /a) 

n n 0 0 


exp|^~ [r^ +r^ - 2rr ^ cos(cp - cp ^ ) 3 | r jdr ^dcp 


= 2 c n I r r 


n 


/2rrrr > 

1 ex _JiO_ 

\ xz . 

' p (x z 


r * +I f +2 *Z 

1 a 1 


i| dr i 


(C. 1) 


The exact evaluation of Eq. (C. 1) is usually impossible, but reason- 
able approximations can be obtained. The Bessel function can be 
asymptotically approximated as 


n 


2 n ! 



- 2^—1 

4(n + 1 ) J 


for x « 1 


V*> e i 


V 



njT TT \ 

2 “ 4 / ' 


for x » 1 


(C.2) 


(C. 3) 


2rrrr j 

For r in the region such that • - ■ — « 1, i. e, , 

A 

will us e Eq. (C.2) in Eq. (C . 1 ) , and 


r « 


X Z 

2tt Tj 


we 
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U(r, cp; z) 


2 C nJ r l [* ' T (~TT^ ) ] exp [i?( r f +r2+ ^a^" r ll 


dr^ . 
(C. 4) 


The point of stationary phase is located at r = -mX Z/a (note that 

X Z 

m can be a negative integer). For r = -mXZ/a , r « 

1 ^ TT r | 

a a _ 1 

= I \ < ~7 , i. e. , the integration in Eq. (C.4) only contributed 

2rr | tn| 6 

to r < a/6 . The area r < a/6 is so small, that it can be neglected. 
Equation (C. 1) can be approximated as 


F/r\ 1/2 / 

U(r,<p;*)=2X J \~) 

n 0 \ 


2 , 2 . 2 X Z 

r + r _ + nr - 

la 1 


cos 


( 2TTrr 1 n\. 

\ TzT 'I) dr i • 


(C.5) 


We can expand the cosine term as the sum of two exponential functions. 
There are two points of stationary phase, one for each exponential term 
of the cosine function, one is r^ = r- nX Z/a , the other is 
r = -r - nXZ/a . Since r^ >0 , not every diffraction order (n) falls 
into the region of integration. We can divide r into a set of zones 
Z(m), given by Z(m) = [mXZ/a, (m + l)XZ/a], m s 0 integer. For a 
given zone Z(m), if m^n then only r + is inside the region of inte- 
gration, if l^(m+l)< -n then only r is inside the region of integration. 
Using the formula for the method of stationary phase, Eq. (C. 5) can be 


written as 
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U(re Z(m),cp ; Z) = 


2 

nS m 


' 2 2 

C expIZrri nr /a - inX Zn / a ] 
n L 



♦ I 

r\<. - {m + 1 ) 


C 

n 


exp 



1/2 

„ . / . , ^ 2/ Z\( , n X Z \ 

2tt mr/a - itt X Zn /a M -1 — — — 1 


(C. 6) 


Changing the index n to -n in the second term of Eq. (C.6), we have 

U (r , V ;Z) = S ) exp(-iTT X Zn 2 /a^) exp(2-ni nr /a) . 

(C. 7) 


For r » XZ/a , then 

U(r,cp; Z) = / C exp(-irrXZn 2 /a 2 ) exp(2ni nr /a) . (C. 8) 

*—*» n 


The result is similar to the Talbot image property of linear 

2 

gratings as discussed in Appendix B* Foi Z - 2Ma /\ f the field is 
the self image of the circular grating* From Eq. (C. 8) we can con- 
clude that each spherical wave propagates independently from the other 

2 , 2 

spherical waves and with relative retarded phase exp(in X Zn / a ) for 
order n, as in the case of plane waves. 



